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ABSTRACT 


The  objective  of  this  research  is  to  create  a  flexible  code  which  is  to  be  used 
in  the  investigation  of  optimuin  (minimum  weight)  shell  designs.  A  shell 
analysis/design  program  (DAPS3)  and  a  general  purpose  numerical  optimization 
program  (ADS)  are  incorporated  into  a  single  code  (THESIS).  This  code  provides 
the  user  great  flexibility  in  changing  the  design  variables  and  constraints  which  model 
the  optimization  problem.  The  optimum  designs  produced  by  this  code  are 
compared  to  DAPS3  optimum  designs  in  order  to  identify  any  improvements  made 
by  the  numerical  optimization  technique. 
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I.  INTRODUCTION 


A.  DEVELOPING  A  MATHEMATICAL  MODEL  FOR  OPTIMIZATION 

Design  optimizatioo  is  a  powerful  tool  available  to  the  design  engineer.  While 
it  is  not  a  new  concept,  onty  recently  has  numerical  design  optimization  become  a 
technique  which  is  routinely  applied  by  practicing  engineers  to  the  myriad  of  design 
tasks.  This  trend  can  be  attributed  to  recent  technological  advances  in  high-speed 
digital  computers,  which  are  necessary  to  solve  the  multitude  of  equations  modeling 
the  optimization  problem. 

One  area  particularly  suited  for  design  optimization  applications  is  structural 
design.  The  design  of  a  structure  requires  a  great  deal  of  forethought  by  the  design 
engineer.  Several  restrictions  (constraints)  on  the  design,  such  as  environmental 
effects,  geometry,  and  material  selection,  must  first  be  determined.  Then,  these 
restrictions  are  applied  to  the  design  process  to  produce  a  specific  output  (objective), 
such  as  maximum  structural  response  or  minimum  cost.  The  goal  of  the  design 
engineer  is  to  design  the  best  possible  structure,  in  terms  of  a  specific  desired  output, 
while  adhering  to  several  expected  constraints  on  the  final  design.  That  is,  the  final 
structure  must  be  the  optimum  design  for  anticipated  conditions. 

In  order  to  achieve  this  optimum  design,  the  design  optimization  process 
requires  that  the  design  task  be  defined  as  a  mathematical  model.  This  model 
consists  of  an  objective  function,  constraint  functions,  and  side  constraint  functions. 
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These  three  functions  are  mathematical  equations  expressed  in  terms  of  design 
variables  and  state  variables.  State  variables  are  fixed  quantities  in  the  mathematical 
equations.  Design  variables  are  those  quantities  which  are  allowed  to  vary  during 
the  optimization  process.  The  objective  function,  which  contains  every  design 
variable,  is  a  single  equation  representing  the  quantity  to  be  optimized  (such  as 
weight,  displacement,  cost,  etc.).  Constraint  functions  are  one  or  more  equations 
which  restrict  the  design  variables;  one  or  more  design  variables  are  contained  in 
each  constraint  function.  Side  constraint  functions,  are  equations  which  limit  the 
upper  and  lower  bound  ranges  of  each  design  variable.  Simply  stated,  the  model 
describes,  mathematically,  what  is  to  be  optimized  (objective)  and  the  limitations 
(constraints  and  side  constraints)  on  design  variables. 

Having  converted  the  design  task  into  a  mathematical  mode],  the  design 
engineer  can  program  the  model.  Using  bis  own  optimization  code,  or,  more  likely, 
available  software,  the  design  engineer  may  now  execute  the  process  of  optimizing 
the  objective  function  and  the  associated  design  variables. 

B.  MODELING  A  CYLINDRICAL  SHELL  FOR  MINIMUM  WEIGHT  DESIGN 

The  design  optimization  task  is  to  design  a  minimum  weight  thin  cylindrical 
shell  subjected  to  external  hydrostatic  pressure.  The  shell  can  be  either  monocoque 
or  ring-stiffened  (by  internal  rectangular  cross-section  rings  placed  at  equidistant 
intervals  within  the  shell).  A  single  isotropic  shell  and  ring  materia]  will  be  used  in 
either  case.  A  more  complex  orthotropic  material  or  a  hybrid  construction  material 
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could  have  been  used  (and  supported  by  the  DAPS3  analysis  mode).  However,  due 
to  time  limitations,  only  an  isotropic  material  will  be  investigated.  This  simpler 
material  should  still  credibly  serve  the  purpose  of  demonstrating  the  application  of 
numerical  optimization  techniques  to  a  minimum  weight  shell  design.  The  optimized 
design  is  also  constrained  by  buckling,  strength,  and  geometiy  considerations.  This 
design  task  is  applicable  to  the  design  of  submerged  cylindrical  structures  under  a 
static  pressure  load.  Examples  might  include  submarine  pressure  hulls  and  torpedo 
cases. 

In  order  to  model  the  above  design  task,  a  shell  analysis  program,  DAPS3 
(Design  and  Analysis  of  Plastic  Shells  version  3),  was  used  to  evaluate  the  buckling 
and  strength  constraint  functions  of  a  thin  shell  under  external  hydrostatic  pressure. 
Development  of  these  constraint  functions  is  discussed  in  Chapter  II  (Shell  Buckling 
and  Strength).  The  optimization  problem  was  evaluated  by  the  general  purpose 
numerical  optimization  program  ADS  (Automated  Design  Synthesis  version  1.10). 
Development  of  the  mathematical  model  for  optimization  is  discussed  in  Chapter  111 
(Optimization  Problem).  Both  programs,  which  are  written  in  the  FORTRAN 
language,  are  incorporated  into  a  single  FORTRAN  code  named  THESIS. 

The  THESIS  code  produces  a  minimum  weight  shell  design,  and  it  indicates  the 
shell  dimension  design  variables  (thickness,  ring  height,  and  ring  width)  associated 
with  that  optimum  design.  Therefore,  the  design  task  is  modeled  and  optimized 
within  the  THESIS  code. 
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C.  THESIS  OBJECTIVES 

DAPS3  uses  an  iterative  technique  to  optiinize  the  shell  design  variables.  The 
variables  are  manipulated  (in  a  "DO"  loop  sense)  until  a  design  is  produced  which 
does  not  buckle  under  the  input  hydrostatic  pressure.  When  all  possible  buckling 
loads  are  greater  than  the  hydrostatic  pressure,  DAPS3  terminates  the  iterations. 
The  DAPS3  generated  optimum  design  is  produced  from  the  last  design  iteration. 
It  was  felt  that  incorporation  of  the  more  sophisticated  numerical  optimization 
technique  should  improve  upon  the  DAPS3  iterative  design  technique. 

In  pursuit  of  the  above  hypothesis,  this  thesis  will  investigate  optimum  designs 
generated  by  the  THESIS  code  over  a  wide  range  of  geometr}-  configurations  (L/OD 
ratios)  and  loads  (external  hydrostatic  pressures).  These  designs  will  then  be 
compared  to  the  optimum  designs  generated  by  DAPS3  (when  used  in  the  "design" 
mode).  Specifically,  the  following  design  elements  will  be  compared;  1)  objective 
(shell  weight),  2)  design  variables  (thickness,  ring  height,  and  ring  width),  and  3) 
constraints  (buckling,  strength,  and  geometry).  Any  design  improvements  made  by 
use  of  the  numerical  optiirlzation  technique  will  be  identified.  Results  are  presented 
in  Chapter  IV.  Conclusions  as  to  the  advantages  and  disadvantages  of  incorporating 
a  numerical  optimizatioti  technique  are  presented  in  Chapter  V. 
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II.  SHELL  BUCKUNG  AND  STRENGTH 


The  purpose  of  this  chapter  is  to  explain  the  formulation  of  the  buckling  and 
strength  (or  stress)  terms  used  by  the  THESIS  code.  A  detailed  derivation  of  the 
buckling  equations  is  provided  by  Renzi  [Refs.  1,2]. 

A.  BUCKLING  EQUATIONS 

Structural  failure  by  buckling  is  associated  with  an  "unstable"  design.  Material, 
geometry,  and  load  factors  contribute  to  this  instability.  The  THESIS  code  designs 
the  geometry  of  a  "stable"  cylindrical  shell  based  on  material  and  load  factors  which 
are  input  by  the  user.  Within  the  THESIS  code,  DAPS3  is  used  as  a  subroutine 
analysis  program  to  calculate  the  buckling  loads  on  either  a  monocoque  or  ring- 
stiffened  shell.  These  buckling  loads  are  then  returned  to  the  numerical  optimization 
program  (ADS),  which  is  incorporated  in  the  THESIS  code,  to  be  considered  in 
buckling  constraints  on  the  optimum  shell  design.  Of  course,  other  constraints 
(strength  and  geometry)  are  also  considered  before  the  THESIS  code  produces  the 
optimum  shell  design. 

The  three  types  of  shell  buckling  loads  considered  are: 

1.  Axisymmetric  Collapse  Pressure  (PC) 

2.  Elastic  Monocoque  Shell  Instability  Pressure  or  Elastic  Interbay  Instability 
Prctsiue  (PIB) 
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3.  Elastic  General  Instability  Pressure  (PGEN) 


Note  that  the  above  buckling  terms  are  expressed  in  terms  of  pressure  rather  than 
in  terms  of  force.  The  specific  unit  of  pressure  used  in  the  THESIS  code  is  (Ib/in^). 
The  use  of  buckling  pressures,  rather  than  buckling  forces,  allows  ea^  comparison 
of  the  buckling  loads  to  the  applied  load,  which  is  the  external  hydrostatic  pressure 
(PSTR),  also  in  (Ib/in^).  The  ccmparison  of  buckling  load  to  applied  load  will 
quickly  reveal  whether  the  shell  will  buckle  under  the  specified  applied  load.  From 
herein,  any  use  of  the  term  "load"  will  mean  pressure  in  units  of  (Ib/io^). 

Since  this  thesis  examines  both  monocoque  and  ring-stiffened  shells,  the  type 
of  shell  also  determines  which  particular  buckling  loads  are  present.  For  the 
monocoque  shell  and  ring-stiffened  shells,  two  buckling  loads  are  present: 
axisymmetric  collapse  pressure  (PC)  and  elastic  shell  monocoque  shell  instability 
pressure  (PIB).  When  the  shell  is  ring-stiffened,  a  third  buckling  load,  elastic  general 
instability  pressure  (PGEN),  is  also  present. 

Note  that  the  buckling  term  PIB  is  associated  with  two  different  names.  If  the 
shell  is  monocoque,  then  PIB  is  caL'ed  the  elastic  shell  monocoque  instability 
pressure.  Otherwise,  the  same  term  is  called  the  elastic  interbay  instability  pressure 
if  the  shell  is  ring-stiffened.  This  distinction  is  discussed  in  the  subsection  concerning 
development  of  the  PIB  term. 
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1.  Axi^'minetric  Collapse  Pressure  (PC) 

This  buckling  load  is  present  in  both  monocoque  and  ring*stiffeoed  shells. 
It  describes  local  arisymmetric  buckling  of  the  shell  in  between  two  ring'stiffeners 
(or  in  between  the  end  frames  for  a  monocoque  shell).  An  example  of  axisymmetric 
buckling  between  two  rings  is  illustrated  in  Figure  1. 


Figure  1.  Local  Axisymmetric  Buckling  of  a  Rlng'Stlffened  Shell  Between 
Adjacent  Rings. 

To  find  an  explicit  expression  for  the  axi^mmetric  collapse  pressure  (PC), 
first  consider  the  differential  element  of  a  thin  cylindrical  shell  shown  in  Figure  2. 
After  establishing  equUibrium  of  forces  and  moments  on  the  element,  it  is  found  that 
the  governing  differential  equation  of  the  isotropic  sheU  is; 
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^  ,  [<gqs(%;|^ « [«rL .  .(i-sTsxi-^) 
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^2  2 
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Figure  2.  DlfTcrentlal  Element  of  an  Isotropic  Thin  Shell 


Subjected  to  External  Hydrostatic  Pressure  (PSTR). 

The  coefficient  of  the  second  order  term,  in  this  fourth  order  ordinary 
differential  equation,  renders  the  solution  w(x)  to  be  nonlinear  with  respect  to  the 
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external  hydrostatic  pressure  (PSTR).  lliis  same  coefficient  also  makes  it  possible 
to  extract  the  axiqrmmetric  collapse  pressure  (PC)  from  the  solution  w(x). 

The  solution  of  the  governing  differential  equation  is  found  to  be: 

H<x)  ■  A  8inh(X,x)  *  B  coA^k^x)  * 

C  »  F  codi<V)  -  *'^**(1~«*) 

However,  due  to  selection  of  the  sxi;:.  -ordinate  x  (or  u)  origin  at  a  point  midway 
between  ring-stiffeners,  two  of  th  arbitrary  constants  are  found  to  be  A  »  C  =  0. 
So,  the  general  solution  reduces  to: 

H<x)  -  B  co8h(A,x)  ♦  F  cosbikjx)  - 

Ed  2 

where: 

X,  ,  X)  are  characteristic  roots 

The  arbitrary  constants  B  and  F  are  solved  using  two  boundary  conditions 
which  require  that: 

1.  The  shell  span  in  between  ring-stiffeners  behaves  like  a  beam  fixed  at  both 
ends. 

2.  Radial  continuity  exists  at  the  ring-stiffener  and  shell  junction. 

These  boundary  conditions  are  applied,  and  it  is  found  that  the  arbitrary  constants 
B  and  F  are  algebraically  quite  complex.  However,  both  constants  possess  the  same 
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denominator  quantity.'  This  denominator  quantity  is  found  to  be  a  function  of  shell 
material  variables,  shell  geometry  variables,  and  the  characteristic  roots.  Since  the 
variables  of  shell  material  and  geometry  are  independent  of  the  PC  term,  the 
axisymmetric  collapse  pressure  (PC)  must  be  extracted  from  the  characteristic  roots. 

The  PC  term  may  be  extracted  from  the  solution  w(x)  by  requiring  that, 
upon  collapse,  the  radial  displacement  of  the  shell  tend  to  infinity.  That  is,  the 
common  denominator  (which  is  a  function  of  the  characteristic  roots)  of  the  arbitrary 
constants  B  and  F  is  equal  to  zero  when  the  external  hydrostatic  pressure  (PSTR) 
exceeds  the  axisymmetric  collapse  pressure  (PC).  Thus: 

The  characteristic  roots  X|  and  Xj  correspond  to\ 

Ml  >1.0 

PC 


Then,  it  is  found  that  the  characteristic  roots  are  pure  imaginary  numbers: 


v/Rf 


PSTH 


-1  + 


PC  N  PC 


PSTR 


♦  1 


y/Rf 


PSTR 


-1  - 


N  PC  N  PC 


PSTR 


+  1 


‘  For  brevity,  the  algebraically  complex  expression  of  this  common  denominator 
is  omitted.  See  [Ref.  2:  p.  A-11]  for  details. 
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When  these  characteristic  roots  are  placed  into  the  common  denominator, 
which  is  then  set  equa]  to  zero,  a  transcendental  equation  results.  This  equation  is 
a  function  of  PC  and  can  be  solved  for  the  same. 

The  DAPS3  subroutine,  in  the  THESIS  code,  calculates  the  PC  term  for 
each  change  in  the  shell  geometry  design  variables  (thickness,  T,  for  monocoque 
shells;  thickness,  T,  ring  height,  HEIGHT,  and  ring  width,  WIDTH,  for  ring>sti£fened 
shells).  Recall  that  the  other  terms  (shell  material  and  geometry  properties  and 
external  hydrostatic  pressure)  used  to  extract  the  PC  term  are  user-specified  and 
remain  constant  throughout  the  process.  The  PC  term  is  then  returned  to  the  ADS 
optimization  program,  incorporated  in  the  THESIS  code,  to  be  evaluated  as  a 
buckling  constraint  term. 

2.  Elastic  Interfoay  Instability  Pressure  (PIB) 

Recall  that  the  PIB  term  is  associated  with  two  different  names: 

•  Elastic  Monocoque  Shell  InstabiUty  Pressure  (for  the  monocoque  shell). 

•  Elastic  Interbay  Instability  Pressure  (for  the  ring-stiffened  shell). 

This  buckling  load  describes  local  non-^rmmetric  collapse  of  a  ring-stiffened  shell  in 
between  any  two  rings.  Similarly,  this  describes  non-^metric  collapse  of  a 
monocoque  shell  in  between  its  two  rigid  bulkheads  (end-rings).  An  example  of 
interbay  buckling  is  illustrated  in  Figure  3.  Regardless  of  shell  type,  from  hereon, 
the  PIB  term  will  be  called  Elastic  Interbay  Instability  Pressure,  which  describes  the 
more  general  case  of  a  ring-stiffened  shell. 
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Figure  3.  Local  Nco*s)r«un€ti1c  Buckling  of  a  Ring>Stifrcned  Shell  Between 
AdlJaecnt  Rings. 

Interbay  buckling  differs  from  the  type  of  buckling  described  by  the  PC 
term  in  subsection  (a)  in  that  inward  and  outward  lobes  are  formed  ahematety 
aro'-Qd  drc«iinfereiice.  That  is,  the  collapse  is  not  axityrometric.  The  number 
of  vir^»f>icrer  iial  waves  formed  by  these  lobes  is  known  as  the  buckling  mode,  which 
is  designated  by  che  variable  (n).  This  buckling  mode  term  will  be  discussed  in  more 
detail  during  development  of  the  equations  necessary  to  find  the  elastic  interbay 
instability  pressure  (PIB). 

Tbr  first  step  in  finding  PIB  is  to  establish  the  equilibrium  equations  for 
a  thin  shell  diffe;cntial  element,  as  shown  is  Figure  4.  The  differential  equflibrium 
equations  -  in  terms  of  forces  (F)),  moments  (M),  and  radial  displacement  (w)  -  in 
the  X,  and  r  directions  are,  respectively: 


Fi^u'.e  4. 


DifTerential  Element  of  an  Isotropic  Thin  Shell 


a*AL. 

►  R  ** 

•  0 

ax* 

d<^ax 

a*//. 

R - - 

dN. 

+  — 2  + 

ax* 

a«i>* 

d^dx 

a*3/. 

ax* 

axd4> 

- i  2  +  AT.  =  + 


R^ 


a4>*  2  ax* 


) 


Next,  a  solution  for  the  displacement  terms  ••  u(axial),  v(circumferential), 
and  w(radial)  ••  is  assumed;* 


*  In  [Ref.  1,  p.2>7],  Renzi  notes  that  this  is  the  same  displacement  pattern  used 
by  Von  Mises  for  isotropic  shelb.  It  is  a  weU  accepted  solution  which  is 
somewhat  conservative. 
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u  •  A  8in(/i4>)siD-^ 
iC 

V  •  B  cos(n^)cos-^ 

R 

w  ~  sin(n<^)cos— 

R 

wfarae; 

P  s  —  ;  a  Unsiqjpocted  shell  length  between  lings 

These  assumed  displacements  account  for  the  sinusoidal  nature  of  the 
circumferential  waves  which  form  upon  collapse.  As  noted  earlier,  the  buckling 
mode  (n)  is  the  number  of  waves  formed. 


Figure  5.  Buckling  Mode  (n)  under  Elastic  Interbay  Instability  Pressure. 


Next,  strain  and  curvature  expressions  are  generated  in  terms  of  the 
assumed  displacements.  Then,  the  forces  (N)  and  moments  (M)  are  expressed  in 


terms  of  the  generated  strain  and  curvature  oppressions.  Subsequently,  the  force  and 
moment  expressions  in  the  established  differential  equilibrium  equations  are  now  in 
terms  of  the  assumed  displacements.  That  is,  the  three  equilibrium  equations  are 
now  in  terms  of  known  shell  material  properties,  known  shell  geometry  properties, 
unknown  arbitrary  constants  A  and  B,  and  the  unknown  term  PIB.  The  following 
system  of  equations  results; 


«4 

k 


Note  that  the  above  a,  terms  are  expressions  in  terms  of  known  shell  geometry 
properties,  known  shell  material  properties,  and  buckling  mode  (n ).  To  insure  that 
the  last  equation  in  the  above  system  is  a  consistent  linear  combination  of  the  other 
two,  the  determinant  of  the  following  augmented  matrix  must  equal  zero; 


det 


1^ 


fl, 


0 


Thus,  the  elastic  interbay  instability  pressure  (PIB)  is  given  by  the 
following  explicit  expression: 
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PIB 


- a, 

«i«5-‘V*4 

The  only  remaining  unknown  in  the  expression  is  the  buckling  mode  (n).  The 
DAPS3  subroutine  in  the  THESIS  code  iterates  the  variable  n  over  an  integer  range 
from  2  to  20.  The  minimum  PIB  calculated  within  this  range  is  the  elastic  interbay 
instability  pressure,  and  the  corresponding  n  is  the  number  of  circumferential  waves 
in  which  the  shell  will  buckle.  The  value  of  PIB  is  returned  to  the  ADS  optimization 
program,  which  is  incorporated  into  the  THESIS  code,  for  evaluation  as  a  buckling 
constraint  term. 

3.  Elastic  General  Instability  Pressure  (PGEN) 

The  third  buckling  load  considered  describes  the  overall  or  general 
instability  of  a  ring*stiffened  shell,  whereas  the  previously  described  PC  and  PIB 
buckling  loads  described  local  instabilities  (in  between  any  two  ring-stiffeners)  for 
both  monocoque  and  ring-stiffened  shells.  Collapse  under  elastic  general  instabih'ty 
pressure  (PGEN)  is  somewhat  similar  to  collapse  under  elastic  interbay  instability 
pressure  (PIB)  in  that  the  shell  does  not  collapse  axitymmewically.  That  is,  the  shell 
collapses  in  inward  and  outward  lobes  around  the  circumference  as  was  seen  in 
Figure  5.  However,  the  buckling  effect  is  not  isolated  between  any  two  ring- 
stiffeners.  Instead,  the  entire  length  of  the  shell  (between  its  rigid  bulkheads) 
collapses  with  a  discernible  buckling  mode  (n),  as  seen  in  Figure  6. 
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In  deriving  the  pre\ious  local  buckling  loads,  PC  and  PIB,  differentia] 
equations  were  used  to  find  a  closed  form  exact  solution.  However,  use  of  a 
differential  equation  method  in  deriving  PGEN  would  be  quite  difficult  without 
making  simplistic,  perhaps  unrealistic,  assumptions.  So,  an  alternative  method  is 
needed  to  obtain  an  e]q>ression  for  PGcN. 

The  approach  used  is  the  Rits  energy  method.  Consider  a  ^tem 
consisting  of  a  ring-stiffened  shell  under  external  hydrostatic  pressure  (PSTR).  Pre¬ 
buckling  displacements  and  post-buckiing  (colliq>se)  displacements  (in  the  shell  and 
in  the  rings)  contribute  strain  energy  to  the  system.  In  addition,  work  is  done  by  the 
applied  load  (PSTR)  to  cause  these  displacements.  Tliua,  the  system  possesses  smne 
total  potential  energy,  which  is  the  sum  of  these  strain  potential  energies  and  the 
work  done  by  external  pressure: 
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where: 


U. 


Total  Potential  Energy 


m  Shell  Extensiooal  Strain  Potential  Energy 
>  Shell  Bolding  Strain  Potential  Energy 
°  Ring  Extewrional  Strain  Potoitial  Energy 
«  Ring  Bending  Strain  Potmtial  Energy 
Wfsnt  *  of  the  External  Prefl8ure(PS770 


Under  the  Ritz  energy  method,  this  total  potential  energy  must  be  a  minimum.  Hiis 
requirement  ultimately  leads  to  a  solution  for  the  elastic  general  instability  pressure 
(PGEN),  which  is  the  critical  PSTR  at  which  the  i^g>stiffened  shell  buckles. 

After  a  long  series  of  derivations,  it  can  be  shown  that  the  potential 
energy,  Uj;  is  a  function  of:  known  shelUring  material  and  geometry  properties;  the 
unknown  elastic  general  instability  pressure  (PGEN);  and  the  variable  displacements 
u,v,  and  w.  The  displacements  due  to  buckling  are  assumed  to  be  of  the  form: 


u  *  i4jCos(n<l>)cos|~j 

V  =  B,8in(/i<|>)siii^-Hj  +  Rj8in(n<»jl-co8|-^| 
w  »  C,cos(/r4>)8in|~j  +  CjCos(/r4>)jl -co6|-^| 
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By  substituting  the  above  displacements  into  the  U,  expression,  the  total 
potential  energy  becomes  a  function  of  the  five  arbitraiy  constants: 

V.  -  ^2. 

Recall  that  the  total  potential  energy  must  be  a  minimum  to  satisfy  equUibrium 
requirements,  thus: 


du,  au,  au,  au,  au, 


The  result  is  five  symmetrical  linear  homogeneous  equations  in  the  unknown 
arbitraiy  constants.  These  equations  can  be  expressed  in  form  (A^]{Xi}  >=  0,  or: 


i4j|  A^2  A^2  i4j4  i4,j 

'll 

^21  -^22  "^23  ^  '*25 

^31  ^32  ^33  ^34  ^35 

<^I 

"^41  -^42  ^43  ^*44  '^45 

^2 

;^51  ^52  ^53  ^*54  ^55, 

Cl 

The  terms  in  the  above  5X5  coefficient  matrix  are  given  by: 

‘  i  ij  • 

where: 

*  expressions  in  terms  of  known  material 
and  geometry  properties  and  the 
budding  mode  va^blefn) 
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Provided  that  the  S  X  5  coefficient  matrix  is  symmetric  and  positive 
definite,  it  has  five  positive  eigenvalues,  (PGEN)i,  and  five  corresponding 
eigenvectors,  {x}^.  The  minimum  eigenvalue  is  the  external  pressure  (PGEN)  at 
which  the  ring'Stiffened  shell  will  buckle. 

However,  the  buckling  mode  (n)  still  remains  a  variable  quantity  to 
consider  when  establishing  the  minimum  eigenvalue.  The  DAPS3  subroutine  iterates 
(n)  in  the  integer  range  fi’om  2  to  12,  and  it  saves  the  minimum  eigenvalue  (of  the 
five  possible  eigenvalues)  for  each  corresponding  mode.  The  array  of  these 
eigenvalues  is  examined,  and  the  absolute  minimum  is  chosen  as  the  elastic  general 
instability  pressure  (PGEN).  This  is  the  value  which  is  returned  to  the  ADS 
optimization  program  to  be  considered  as  a  buckling  constraint  term  for  a  ring* 
stiffened  shell. 

B.  STRENGTH  (STRESS)  EQUATIONS 

Another  type  of  failure  to  consider  is  that  of  strength  limitation.  By  this,  it  is 
meant  that  the  shell  fails  when  resultant  Von  Mises  stresses,  at  specific  locations, 
exceed  the  yield  strength  of  the  sbelVring  material.  Since  material  properties  (such 
as  elastic  modulus,  Poisson’s  ratio,  yield  strength,  etc.)  are  input  by  the  THESIS  code 
user  and  remain  fixed,  some  constraints  must  be  placed  on  the  optimum  design  to 
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account  for  material  limitations.  These  limitations  are  the  strength  (or  stress) 
constraints.’ 

Since  these  strength  constraints  compare  q>ecific  stresses  to  the  materia]  yield 
strength,  it  is  important  to  select  appropriate  locations  at  which  stresses  are  to  be 
analyzed.  There  are  four  shell  locations  chosen  as  representative  points  most  likely 
to  fail  by  yield: 

1.  The  middle  fiber  of  the  shell  at  mid-bay  (half  way  between  any  two  ring- 
stiffeners  or  the  two  rigid  bulkheads  of  a  monocoque  shell).  This  stress  is 
designated  by  the  variable  STR. 

2.  The  outer  fiber  of  the  shell  at  the  shell  and  ring  junction.  This  stress  is 
designated  by  the  variable  SIGVOF. 

3.  The  inner  fiber  of  the  shell  at  the  shell  and  ring  junction.  This  stress  is 
designated  by  the  variable  SIGVIF. 

4.  The  center  fiber  of  a  rectangular  cross-section  ring.  This  stress  is  designated 
by  the  variable  SRF. 

The  above  locations,  shown  in  a  cross-section  view  of  a  ring-stiffened  circular 
cylindrical  shell,  are  illustrated  in  Figure  7. 

'fbe  following  subsections  will  briefly  explain  the  development  of  the  equations 
necessary  to  calculate  the  above  Von  Mises  stresses.  The  effects  of  stress 
concentration  are  not  considered  in  the  calculation  of  these  stresses.  Note  that  both 
material  yield  strength  and  calculated  stresses  are  in  units  of  (Ib/in^). 


’  Although  the  terms  "strength  constraint"  and  "stress  constraint"  are 
interchangeable,  the  single  term  "strength  constraint"  will  be  used  throughout 
this  discussion  for  the  sake  of  consistency. 
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ngure  7.  Location!  of  Strtn  Value!  STR,  SIGVOF»  SIGVIF,  and  SRF 
Shown  in  a  Shell  Cro!!>section  View. 

1.  Middle  Fiber  Shell  Stre!!  at  Mid-bay  (STR) 

The  effect!  of  the  radial  itrea!  component  are  ignored  in  the  calculation 
of  this  Von  Mises  streu.  In  fact,  the  radial  component  ia  ignored  in  the  calculation 
of  the  Von  Miaes  atresa  at  all  four  location!  aince  the  ahell  ia  aaaumed  to  be  thin. 
(A  aide  conatraint  ia  included  in  the  THESIS  code  to  inaure  that  the  thickoeaa,  T,  ia 
leaa  than  or  equal  to  20  percent  of  the  mean  radiua,  R.  Tbia  inaurea  that  the 
optimum  design  has  negligible  radial  stress  components.) 
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Calculation  of  the  middle  fiber  stress  value,  STR,  begins  with  a  definition 
of  the  Von  Mises  stress  at  the  inner  and  outer  fibers  (at  mid-bay).  The  Von  Mises 
stress  of  these  fibers  is: 


srfkete: 

■  Von  Mises  stre«  at  mid-bay 
Oj^  •  Axial  stress  at  mid-bay 

■  Circumferential  stress  at  mid-bay 

The  subscript  "i"  and  the  superscript  "o”  correspond,  respectively,  to  inner  and  outer 
fiber  locations. 

Next,  the  mid-bay  Von  Mises  stresses  at  inner  and  outer  fiber  locations  are 
each  divided  by  the  material  yield  strength  to  obtain  two  percentage  values.  The 
middle  fiber  stress  value,  STR,  is  an  average  of  these  two  percentages; 


adiete: 

o,  ■  Material  yield  strength 


The  superscripts  "i"  and  "o"  correspond,  respectively,  to  inner  and  outer  fiber 
locations. 

The  stress  value  STR  is  calculated  by  the  DAPS3  subroutine  in  a  much 
more  complex  scheme  than  described  above;  the  process  of  calculating  the  hoop  and 


23 


axial  stresses  is  quite  involved/  The  STR  value  is  returned  to  the  ADS  optimization 
program  to  be  considered  as  the  first  of  four  strength  constraints  on  the  optimum 
design. 

2.  Outer  Fiber  Shell  Stress  at  a  Ring  Junction  (SIGVOF) 

This  Von  Mises  stress  is  calculated  by  a  method  similar  to  that  of  the  STR 
value  calculation.  The  onty  difierence  being  that  outer  fiber  shell  stress  components 
(at  the  shell/ring  junction)  are  used,  instead  of  the  outer  and  inner  fiber  shell  stress 
components  (at  mid-bay).  The  Von  Mises  outer  fiber  shell  stress  at  the  junction  is; 

-  /K/f  -  K/K”)  *  K/f 

stdiexe: 

Oyyj  ■  Von  Mises  stress  at  shelVting  junction 

Ojij  ■  Axial  stress  at  sheUMng  Junction 

o^j  •  Circumferential  stress  at  shell/ring  junction 

The  "o"  superscript  indicates  an  outer  fiber  location. 

The  strength  constraint  requires  that  this  Von  Mises  stress  be  compared 
to  the  material  yield  strength.  Thus,  the  value  of  SIGVOF  is* 


*  See  Chapter  1  of  (Ref.  1]  for  details. 
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^  # 

SICVOF  - 

This  is  the  value  retuined  to  the  optimizatioii  program  for  consideration  as  the 
second  of  four  strength  constraints. 

3.  Inner  Fiber  Shell  Stress  at  a  Ring  Junction  (SIGVIF) 

The  development  of  this  value  is  identical  to  that  of  the  SIGVOF  value 
except  that  inner  fiber  shell  stress  components  are  used.  The  Von  Mises  inner  fiber 
shell  stress  is  defined  as: 

■  'IM  -  M-’xA  *  KT 

where: 

Oy^j  ■  Von  Misea  stress  at  shell/ring  juncti<m 
a^j  •  Axial  stress  at  sheil/ring  junction 

■  Circumferential  stress  at  shell/ting  junction 


Then,  this  Von  Mises  stress  is  compared  to  the  material  yield  strength: 


SIOVIF  - 


This  is  the  value  returned  to  the  optimization  program  for  consideration  as  the  third 
of  four  strength  constraints. 
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4.  Center  Fiber  Ring  Stress  (SRF) 

The  last  of  the  four  strength  constraints  is  the  Von  Mises  stress  at  the 
centroid  of  a  ring-stitfener.  Note  that  this  is  not  a  shell  stress,  but  rather  a  stress 
inside  the  internally  located  rectangular  cross*8ection  ring.  Since  the  ring'Stiffener 
does  not  bear  an  axial  load,  the  Von  Mises  stress  inside  the  ring  reduces  to  the 
circumferential  stress,  which  is  defined  as: 

o  ‘  .  o  ‘ 

^  ♦  iWIDTH)(J)\KR] 

where: 

Q*  *  Total  load  on  ring  per  unit  circumference  length 

RR  •  Radius  from  shell  centaline  to  ting  centroid 
HEIGHT  «  Height  of  rectangular  ring -stiffener 
WIDTH  •  Width  of  rectangular  ring-stiffener 

The  V  superscript  indicates  a  center  fiber  location,  which  is  the  centroid  of  the  ring. 
Perhaps,  the  above  shell  geometry  terms  are  best  iUustrated  by  Figure  8. 
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Figure  8. 


Illustration  of  Various  Shell  Geometiy  Terms. 


III.  OPTIMIZATION  PROBLEM 


Now  that  the  buckling  and  strength  tenns  have  been  defined,  a  mathematical 
model  of  the  optimization  problem  can  be  defined.  The  model  consists  of  an 
objective  function  (F),  constraint  functions  (^),  and  side  constraint  functions.  These 
functions  will  be  defined  in  subsequent  subsections  of  this  chapter.  A  description  of 
the  Sequential  Linear  Programming  (SLP)  technique  used  to  optimize  the  model 
follows  those  subsection*^ 

A.  OBJECri\'E  FUNCTION 

Recall  that  the  objective  of  the  optimization  is  to  design  a  minimum  weight 
circular  cylindrical  shell.  The  objective  function  (F),  then,  is  a  mathematical 
expression  of  the  shell  weight.  For  the  general  case  of  a  ring>stiffened  shell,  the 
objective  function  is; 

whoe: 

p  ■  Specific  weight  [ibf/ft^j  of  the  sheUfakin) 

Py  ■  Specific  weight  [ibf/ft’j  of  the  ringfs) 

Figure  8  best  illustrates  the  above  variable  terms.  The  first  quantity  represents  the 
shell(skin)  weight,  while  the  second  quantity  represents  the  ring(s)  weight.  Note  that 
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the  weight  is  expressed  on  a  (Ib/iFt)  basis.  The  above  expression  is  also  valid  for  the 
more  pacific  case  of  a  monocoque  shell.  In  that  case,  the  AR  (ring  cross-section 
area)  term  would  just  be  zero. 

The  objective  function  contains  gh  the  design  variables  of  the  optimization 
problem.  Recall  that  the  single  design  variable  for  the  monocoque  shell  is: 

1.  Thickness  (T) 

In  contrast,  the  ring-stiffened  shell  has  three  design  variables,  which  are: 

1.  Thickness  (T) 

2.  Ring  height  (HEIGHT) 

3.  Ring  width  (WIDTH) 

The  design  variable  T  is  explicitly  expressed  in  the  objective  function.  The  objective 
function  (F)  also  contains  indirect  expressions  which  incorporate  the  design  variables. 
These  expressions  are: 

R  *  RO  -  - 
2 

YR  ~  T  * 

2 

AR  -  (HEIGHTKWIDTH) 

Now  that  the  objective  has  been  defined  in  terms  of  the  design  variable(s),  the  first 
step  in  defining  the  optimization  model  is  complete.  Next,  constraints  on  the  design 
variables  will  be  defined. 
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B.  CONSTRAINT  FUNCTIONS 

The  constraints  on  the  design  variables  are  divided  into  three  groups:  1) 
buckling,  2)  strength,  3)  geometry.  As  noted  in  Chapter  II,  the  number  and  type  of 
constraints  placed  on  the  design  variables  depends  upon  the  particular  shell  being 
designed:  monocoque  or  ring-stiffened.  The  next  two  subsections  define  the  specific 
constraints  placed  on  monocoque  and  ring-stiffened  shell  designs. 

1.  Monocoque  Shell  Constraints 

The  monocoque  shell  is  subject  to  two  buckling  loads  (PC,PIB);  four  stress 
loads  (STR,SIGVOF,SIGVIF,SRF);  and  no  geometry  constraints.^  A'xordingly,  the 
six  constraints  on  the  monocoque  shell  design  are: 

»  l-(PCIPSm  s  0 
*  l-iPIBfPSTR)  S  0 
^3  -  (STTVIOO)-!  s  0 
g^  =  (SIGVOF/m)  s  0 
g^  =  (siGViFimyi  i  0 
g^  =  (SRFimyi  s  0 

In  the  above  expressions,  gj(j=l,2)  represent  the  buckling  constraints  and 
gjJ=3,4,5,6)  represent  the  strength  constraints  on  the  monocoque  shell  design. 


^  Actually,  the  monocoque  shell  is  subject  to  a  type  of  geometry  constraint  which 
limits  the  range  of  the  design  variable,  thickness  (T),  between  an  upper  and 
lower  bound.  However,  it  will  be  seen  that  this  is  actually  a  side  constraint  on 
the  design  variable. 
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2.  Ring'StitTened  Shell  Constraints 

The  ring-stiffened  she!J  is  subject  to  three  buckling  loads  (PC,  PIB, 
PGEN),  four  stress  loads  (STR,  SIGVOF,  SIGVIF,  SRF),  and  one  geometry 
constraint  which  limits  (be  ring  aspect  ratio  (HEIGHT/WIDTH)  to  be  less  than  4:1.‘ 
This  geometry  constraint  is  necessary  to  ensure  that  the  resultant  ring 
design  ic  not  "flimsy",  that  is,  unstable.  The  predominate  mode  of  instability  for 
rectangular  ring-stiffcners  is  the  tendency  to  rotate  about  the  ring-stiffener  base  [Ref. 
3:  pp.  23  -24].  This  rotation  causes  the  ring-stiffener  to  deform  in  a  sinusoidal  pattern 
which  crosses  in  and  out  of  the  ring’s  vertical  plane  of  symmetry.  Since  DAPS3  does 
not  check  for  this  instability  when  designing  a  ring-stiffened  shell,  this  particular 
constraint  keeps  the  ring  aspect  ratio  low  enough  that  the  ring  instability  mode  is 
adequately  addressed  (by  the  THESIS  produced  shell  design).^  Inclusion  of  this 
constraint  may  result  in  a  weight  penalty  since  a  "flimsy"  ring  may  weigh  less  than  a 
short  and  thick,  but  inherently  stable,  ring. 

The  resulting  eight  constraints  on  the  ring-stiffened  shell  design  are  given 
by  the  following  gj  expressions: 


^  As  explained  in  the  previous  footnote,  side  constraints  on  the  design  variables 
(T,  HEIGHT,  WIDTH)  also  apply  for  the  ring-stiffened  shell. 

^  DAFS3  does  attempt  to  prevent  a  "flim^'  ring.  If  intermediate  designs  have 
an  aspect  ratio  greater  than  4:1  and  a  ring  width  less  than  0.125  inch,  a 
minimum  ring  width  of  0.125  inch  is  then  forced.  However,  this  does  not 
necessarily  prevent  the  aspect  ratio  from  increasing  to  an  unreasonable  number 
(>>4:1). 
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g^  -  i-{pc/psm  s  0 
»  l-(PIBiPSm  i  0 
»  l-iPGEN/PSm  i  0 
g^  -  (iT^lOC,  I  s  0 
g,  -  iSIGVOFim)  £  0 
g^  »  iSIGVIFIl(JO)-l  £  0 
g^  -  (5W7100)-1  <  0 

-  [(ff£/(?^»7/)7W)/4)-l  s  0 

In  the  above  expressions,  ^(j- 1,2,3)  represent  the  buckling  constraints,  i 

represent  the  strength  constraints,  and  gj(j«8)  represents  the  geometry  consc^int. 


C.  SIDE  CONSTRAINTS 

Side  constraints  directly  impose  bouu6  on  the  design  variables.  The  ADS 
optimization  program  incorporated  into  the  THESIS  code  requires  that  the  user 
specify  lower  and  upper  bounds  on  all  declared  design  variable's.  The  optimization 
search  will  then  only  sxx:ur  in  that  portion  of  the  design  space  delineated  by  the  side 
constrai  itr . 

The  moDocoque  shell  has  only  vce  design  variable,  T,  and  it  is  directly 
constrained  by  the  following  bounds; 

0.10  IN  s  r  s  i.818  IN 

The  lower  bound  was  chosen  as  an  aihitrarily  low  number  to  permit  a  broad  design 
space,  knowing  that  the  lower  bound  would  not  be  approached  at  the  higher  external 
hydrostatic  pressure  shell  loads.  The  upper  bound  was  chosen  such  that  the  shell 
thickness  (T)  would  not  exceed  the  ratio; 
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^^lere; 


^  iOJ20 
R 


R  •  —  {Shell  Outside  Diameter)  -  — 
2  2 


T 

•  RO  -  - 
2 


Since  the  shell  outside  diameter  is  fixed  input  (OD  =  20  IN),  the  outside  radius 
remained  fixed  (RO  ^  10  in).  Thus,  the  value  of  T  «  1.818  IN  was  calculated  such 
that  the  shell  thickness  never  exceeded  20%  of  the  mean  radius.  Recalling  that  the 
objective  is  to  design  a  thin  circular  cylindrical  shell,  a  T/R  value  less  than  20% 
assures  that  the  shell  remains  "thin".  That  is,  radial  stress  is  much  smaller  than  axial 
or  circumferential  stresses  and  can  be  ignored. 

The  design  variables  of  a  ring-stiffened  shell  are  assigned  the  following  side 
constraints: 


0.10  IN  s  r  s  1.818  IN 
0.05  IN  s  HEIGHT  i  10.0  IN 
0.05  IN  s  WIDTH  s  5.0  IN 


Again,  the  lower  bounds  are  chosen  to  be  arbitrarily  low.  The  upper  bounds  are 
chosen  such  that:  HEIGHT  does  not  exceed  the  physical  dimension  of  the  shell 
outside  diameter,  and  WIDTH  does  not  become  so  wide  that  the  resulting  HEIGHT 
is  negligible  (;n  other  words,  the  ring  would  be  virtually  flat,  and  a  monocoque  shell 
would  probably  be  more  appropriate). 
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D.  SUMMARY 


The  mathematical  models  of  the  monocoque  and  the  ring-stiffened  shell  are 
summarized  below.  These  are  the  models  programmed  into  the  THESIS  code  to 
solve  the  optimization  problem. 

1.  Monocoque  Shell 


Minimize; 


F 


24a 

L 


(R)(7XI)(P)  ♦  {RO-YRiiA 


Subject  to: 


g,  -  i-iPQPSm  i  0 
gj  -  \~(JPIBIPSm  s  0 
gj  «  (STTVIOO)-!  s  0 
g<  =  (S/GKOf/100)  s  0 
g,  -  (5/GWF/100)-l  s  0 
gj  -  (SRF/I00)-1  s  0 


Side  constraints; 


0.10  IN  s  r  s  1.818  IN 
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2.  Rlng-SUirened  Shell 


Minimize: 

(WXiKrt « (iK>-no(W-^-ij(p,)] 


Subject  to: 


g^  -  l-(PCIPSTK^  i.  0 
gj  «  \-iPIBIPSm  i  0 
g,  -  \~iPGENIPSm  4  0 

-  (STTVIOO)-!  4  0 
gj  -  {SIGVOFIIQO)  4  0 

«  (S/GWF/100)-1  4  0 

-  (S«7100)-l  4  0 

g,  -  [(H£/t7^ira)7W)/4]-l  s  0 


Side  constraints: 

0.10  IN  4  r  4  1.S18  IN 
0.05  IN  4  HEIGHT  4  10.0  IN 
0.05  IN  4  WIDTH  4  5.0  IN 


E.  OPTIMIZATION  TECHNIQUE 

The  mathematical  models  are  now  structured  for  application  of  an  optimization 
technique.  The  general  purpose  numerical  optimization  program,  Automated  Design 
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Synthesis  (ADS),  is  used  in  conjunction  with  the  analysis  mode  of  the  DAPS3  to 
produce  optimum  shell  designs.  An  excellent  reference  on  numerical  optimization 
techniques  is  found  in  [Ref.  4]. 

The  procedure  to  obtain  the  optimum  shell  design  involves  a  mathematical 
search  frcnn  an  initial  design,  X”,  to  the  optimum  dengn,  X*.  A  gradient  method  is 
used  to  search  for  X*  within  the  design  space  q)ecified  by  the  side  constraints. 
Search  limitations  within  this  design  sp»ee  are  imposed  by  the  constraint  functions. 
The  objective  function,  F,  when  evaluated  at  the  optimum  design,  X*,  is  the 
minimized  objective:  shell  weight. 

The  above  optimization  procedure  is  described  in  more  detail  by  the  following 
general  steps  [Ref.  5;  p.  IS]: 

1.  User  inputs  an  initial  design,  X^  for  the  iteration;  q  »  0. 

2.  Iteration  number,  q,  is  updated  to  q  <■  q  +  1. 

3.  Objective  function,  F(X),  and  the  constraint  functions,  gj(X),  are  evaluated  at 
the  design  X*^'. 

4.  Active  constraints,  J,  are  identiHed.  These  are  the  g^  constraints  which  are 
very  close  to  zero,  within  a  specified  tolerance,  say  0.003. 

5.  The  objective  function  gradient,  V  F(X'’'*)  is  evaluated,  and  critical  constraint 
function  gradients,  V  gi(X^*)  are  evaluated  for  j  e  J. 

6.  Search  direction.  S'*,  is  determined  from  the  evaluations  in  step  5. 

7.  A  one-dimensional  search  is  performed  to  find  the  step  length  a. 

8.  The  design,  X**"*,  is  updated  by:  X’  ■  X’’’  +  o'S". 
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9.  Convergence  criteria  are  examined.  If  design  has  converged,  then  X*  has 
been  found.  Otherwise,  the  procedure  is  reiterated  from  step  2. 

Execution  of  the  above  general  procedure  may  be  accomplished  by  a  variety  of 
methods,  especially  at  steps  3,  6,  and  7.  Step  3  is  termed  the  "strategy"  (ISTRAT); 
step  6  is  termed  the  "optimizer”(IOPT);  and  step  7  is  termed  the  "one-dimensional 
search"(10NED).  These  are  the  three  basic  levels  at  which  the  ADS  numerical 
optimization  program  can  be  controlled  to  produce  an  optimum  design. 

ADS  offers  a  great  number  of  combinations  of  these  basic  levels.  In  fact, 
version  1.10  offers  nine  strategies,  five  optimizers,  and  eight  one-dimensional  search 
methods  [Ref.  6].  Selecting  an  appropriate  combination  of  ISTRAT,  lOPT,  and 
lONED  is  a  matter  of  experience,  since  not  all  possible  combinations  are  suitable 
to  particular  kinds  of  optimization  problems. 

The  specific  combination  chosen  for  this  study  is: 

•  ISTRAT  <=  6;  Sequential  Linear  Programming 

•  lOPT  =  5:  Modified  Method  of  Feasible  Directions 

•  lONED  e  6:  Golden  Section  +  Polynomial  1-D  Search 

The  Sequential  Linear  Programming  (SLP)  method  is  an  optimization  strategy  that 
approximates  a  nonlinear  problem  as  a  linear  problem.  This  is  done  by  creating  a 
first  order  Taylor  Series  approximation  of  the  objective  and  constraint  functions. 
Now  that  the  problem  is  linear,  an  evaluation  of  the  objective  and  constraint 
functions  are  easily  and  inexpensively  calculated,  in  terms  of  computer  time.  Then, 
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the  Modified  Method  of  Feasible  Directions  is  used  to  determine  a  search  direction 
S^.  This  is  also  easily  obtained  since  the  gradients  of  the  objective  and  constraint 
functions  are  directly  available  from  the  Taylor  Series  expression.  Next,  the  Golden 
Section  +  Polynomial  one-dimensional  search  method  is  used  to  find  the  step  length, 
a.  This  step  length  is  how  far  the  search  direction  vector,  S^,  vrill  traverse  from 
design  point  An  c^timum  design,  X^,  to  this  linear  approximation  is  then 
obtained.  The  original  nonlinear  problem  is  again  linearized  about  this  design.  The 
process  is  repeated  until  an  optimum  solution  is  obtained,  as  determined  by  specific 
convergence  criteria. 

This  combination  using  the  SLP  method  is  a  very  fast  optimization  technique 
since  the  original  nonlinear  problem  is  approximated  as  a  linear  problem.  Speed 
represents  the  major  advantage  in  using  this  method.  Possible  disadvantages  include: 
1)  the  optimum  design  may  be  infeasible  and  2)  the  method  may  perform  poorly  for 
an  underconstrained  problem  (Ref.  4:  pp.  156-157].  The  first  disadvantage  means 
that  X*,  though  very  close  to  the  true  optimum,  may  have  a  small  degree  of 
constraint  violation.  The  second  disadvantage  means  that  if  there  are  fewer  active 
constraints  at  the  optimum  than  there  are  design  variables,  the  method  may  not 
rapidly  or  accurately  converge  to  the  true  optimum.  ADS  deals  with  this  difficulty 
by  imposing  finite  move  limits  on  the  linear  approximation  so  that  the  true  optimum 
will  be  approached  within  the  tolerance  of  the  move  limits.  These  possible 
disadvantages  are  not  considered  severe  enough  to  preclude  use  of  the  .SLP  method 
for  this  optimization  problem.  In  fact,  all  methods  within  the  ADS  program  have 
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certain  advantages  and  disadvantages  which  must  be  accepted  when  that  method  is 
applied  to  a  particular  optimization  problem. 

The  optimization  combination  chosen  for  this  study  was  compared  against  other 
direct  and  indirect  methods.*  It  was  found  that  the  array  of  optimum  designs 
produced  by  the  SLP  method  were  almost  identical  to  the  other  strategies.  The 
q>eed  advantage  offered  by  the  SLP  method  was  the  principal  reason  that  this 
particular  optimization  combination  was  chosen.  The  next  chapter  presents  the 
results  of  applying  SLP  to  the  mathematical  models  defined  in  this  chapter. 


*  An  "indirect"  method  involves  converting  a  constrained  problem  into  an 
unconstrained  problem;  the  constraints  are  accounted  for  by  a  penalty  function 
method,  and  the  objective  function  is  converted  into  a  pseudO'objective 
function  which  includes  the  penalty  function.  By  contrast,  a  "direct"  method 
involves  dealing  with  the  objective  and  constraint  functions  individually  during 
the  optimization  process.[Ref.  4:  pp.S4'58] 


39 


IV.  RESULTS 


A.  MONOCOQUE  SHELL 

Tbe  elements  of  an  q|)timum  design  are  derived  from  the  optimization  model. 
Accordingly,  the  results  are  presented  in  three  parts: 

1.  Objective:  DAPS3  vs.  THESIS 

2.  Design  Variables:  DAPS3  vs.  THESIS 

3.  Constraints:  DAPS3  vs.  THESIS 

Tbe  first  part  will  compare  the  optimum  shell  weights  achieved  by  each  program 
(DAPS3  and  THESIS).  The  second  part  then  compares  the  optimum  shell  thickness’ 
(T)  which  correspond  to  these  optimum  shell  weights.  Tbe  last  part  will  compare 
constraint  activity  (buckling  and  strength  violations)  of  tbe  two  optimum  designs. 

1.  Objective:  DAPS3  vs.  THESIS 

A  graphical  comparison  of  tbe  optimum  weights  achieved  by  DAPS3  and 
THESIS  is  shown  in  Figures  9  and  10.  The  graphs  are  set  up  such  that  tbe  optimum 
monocoque  shell  weights  are  shown  over  a  range  of  external  hydrostatic  pressure 
(PSTR)  loads  and  over  a  range  of  shell  geometries  (L/OD  ratios,  which  were  varied 
from  9:1  to  1:1).  It  can  be  seen  that  tbe  weight  curves  (of  each  L/OD  ratio)  is 
essentially  the  same  for  both  programs.  The  one  noticeable  difierence  is  seen  in  tbe 
1:1  (L/OD  ratio)  weight  curves.  The  THESIS  objective  (weight)  is  somewhat  higher 
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WEIGHT  [LB/FT]  J®  WEIGHT  [LB/FT] 


DAPS3  Optimum  Shell  Weight 


THESIS  OPTIMUM  TOIGHT 


PSTR  [PSI] 


Figure  10.  THESIS  Optimum  Shell  Weight 
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than  the  DAPS3  objective  after  PSTR  =  1000  psi.  This  will  be  accounted  for  in  part 
three  of  these  results.  It  will  be  seen  that  the  increase  in  weight  was  necessary  due 
to  an  increase  in  thickness  to  compensate  for  a  constraint  violation. 


At  this  point,  the  results  of  an  optimum  shell  weight  comparison  do  not 
indicate  the  advantage  of  using  a  numerical  optimization  technique.  The  next 
optimum  design  element  to  be  examined  is  the  single  design  variable  of  the 
monocoque  shell:  shell  thickness  (T). 

2.  Design  Variables:  DAPS3  v&  THESIS 

The  design  variable,  shell  thickness  (T),  is  a  monotonically  increasing 
variable.  That  is,  as  the  external  hydrostatic  pressure  (PSTR)  increases,  the 
thickness  also  increases.  This  results  in  the  shell  weight  being  directly  proportional 
to  the  shell  thickness,  as  would  be  expected.  Accordingly,  the  optimum  thickness 
curves  look  like  the  optimum  weight  curv  2S.  This  can  be  seen  in  Figures  11  and  12. 

Since  the  monocoque  shell  only  has  one  design  variable,  the  design  space 
is  quite  simple.  There  is  only  one  design  variable  which  corresponds  to  the  global 
optimum  in  this  design  space.  Even  the  simple  iterative  technique  employed  by 
DAPS3  will  converge  quite  close  to  the  unique  design  variable  which  corresponds  to 
the  global  optimum  (minimum  weight).  Thus,  the  advantage  of  using  a  numerical 
optimization  technique  to  find  the  optimum  design  variable  is  not  realized.  The  last 
part  of  these  results  will  indicate  the  advantage  of  using  numerical  optimization  in 
even  a  simple  one  variable  design  space. 
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DAPS3  OPTIMUM  THICKNESS 

MONOCOQUE  SHELL 


Figure  11.  DAPS3  Optimum  Sheil  Thickness  (T) 


THESIS  OPTIMUM  THICKNESS 


MONOCOQUE  SHELL  i/od 
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2000  3000  4000 
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0000 

PSTR  [PSI] 

Figure  12.  THESIS  Optimum  Shell  lliickness  (T) 
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3.  Constraints:  DAPS3  vs.  THESIS 

The  six  constraints  imposed  on  a  monocoque  shell  design  are  shown  in 
Figures  13  through  18.  The  graphs  are  set  up  such  that  a  violated  constraint  is  a 
positive  gj(j  s  1,6)  value;  a  satisfied  constraint  is  represented  by  a  negative  value; 
and  an  active  constraint  is  a  zero  value.  If  the  constraint  is  violated,  this  is  an 
indication  of  an  infeasible  design.  A  feasible  design  is  indicated  by  an  active 
constraint  or  a  satisfied  constraint. 

The  purpose  of  showing  graphs  of  each  constraint  is  to  identify,  at  a 
glance,  any  violations  which  the  "optimum"  design  may  have.  The  optimum  designs 
produced  by  DAPS3  and  THESIS  may  be  a  minimum  weight  shell  design  for  the 
input  load  (PSTR)  and  geometry  (L/OD).  However,  the  code  user  must  be  able  to 
verify  that  the  produced  optimum,  in  fact,  satisfies  the  original  constraints  imposed 
upon  the  design.  If  there  are  violations,  the  design  produced  is  not  a  valid  design. 
That  is,  the  design  is  infeasible:  the  objective  may  have  been  minimized,  but  at  the 
expense  of  violated  constraints. 

In  comparing  DAPS3  and  THESIS  objective  values,  one  program’s  design 
may  weigh  less  than  the  other  (when  load  r'.nd  geometry  are  the  same  for  each 
program).  Before  accepting  the  lesser  weight  dcsivo  as  the  tme  r^trmum,  constrakt; 
activity  must  be  verified  as  outlined  above.  An  example  of  this  was  seen  in 
monocoque  shell  designs  of  1:1  L/OD  ratios.  DAFS3  appears  to  produce  lesser 
weight  designs  above  PSTR  -  lOOG  psi  than  THESIS,  however,  these  DAPS3  designs 
violate  the  strength  (stress)  constraints  g)(STR)  and  g5(SIGVIF). 
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DAPS3  G(1 )  CONSTRAINT 

MONOCOQUE  SHELL  u/tio 


Figure  13.  DAPS3  Constraint:  Symmetric  Local  Buckling  (PC) 


THESIS  G(l)  CONSTRAINT 

MONOCOQUE  SHELL  i/op 


Figure  14.  THESIS  Constraint’  Symmetric  Local  Buckling  (PC) 
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Figure  15.  DAPS3  Constraint:  Non^Symmetric  Loi:aI  Buckling  (PIB) 


THESIS  G(2)  CONSTRAINT 

MONOCOQUE  SHELL  i/c 


THESIS  Constraint:  Non-Symmetrlc  Local  Buckling  (PIB) 


DAPS3  G(3)  CONSTRAINT 

MONOCOQUE  SHELL 


Figure  17.  DAPS3  Constraint:  Stress  (STR) 


THESIS  G(3)  CONSTRAINT 

MONOCOQUE  SHELL  i/ot> 


0  1000  2000  3000  4000  5000  6000 


PSTR  [PSI] 

Figure  18.  THESIS  Constraint:  Stress  (STR) 
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DAPS3  G(4)  CONSTRAINT 

MONOCOQUE  SHELL 


PSTR  [PSI] 

Figure  19.  DAPS3  Constraint:  Stress  (SIGVOF) 


THESIS  G(4)  CONSTRAINT 

MONOCOQUE  SHELL  VW) 


Figure  20.  THESIS  Constraint  Stress  (SIGVOF) 
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DAPS3  G(5)  CONSTRAINT 

MONOCOQUE  SHELL  VW 


Figure  21.  OAPS3  Conitraint:  Stress  <'S  ' 


THESIS  G(5)  CONSTRAINT 

MONOCOQUE  SHELL  VOd 


Figure  22.  THESIS  Constraint:  Stress  (SIGVIF) 
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DAPS3  G(6)  CONSTRAINT 

MONOCOQUE  SHELL  i/n 


Figure  23.  DAPS3  Constraint  Stress  (SRF) 


The  DAPS3  stren^  constraint  violations  are  seen  in  Figures  17  and  21. 
This  is  especially  noticeable  in  Figure  21  T1  .  nh  constraint  (SIGVIF) 

remains  violated  at  all  pressures  above  JOD)  and  above  4000  psi 

(for  3:1  IVOD).  Thus,  all  1:1  and  3:1  DAPS3  designs  above  the 

aforementioned  pressures  are  infeasible  due  to  matt  rial  failure  by  yield.  By  contrast, 
the  THESIS  designs,  which  include  strength  constraints,  designs  monocoque  shells 
which  satisfy  all  strength  constraints.  This  can  be  seen  in  Figures  18,  20,  22,  and  24. 

The  DAPS3  iterative  technique  is  only  constrained  by  buckling 
considerations.  As  seen  in  Figures  13  and  IS,  the  DAPS3  designs  do  not  violate 
either  monocoque  shell  buckling  constraint  throughout  the  range  of  applied  loads. 
DAPS3  designs  a  minimum  weight  shell  with  the  consideration  being  that  a  shell 
which  just  exceeds  buckling  loads  will  be  the  absolute  minimum  design  required  to 
withstand  the  external  hydrostatic  pressure.  All  monocoque  designs  thicker  than  this 
absolute  minimum  will  also  not  buckle.  These  designs  will  also  withstand  greater 
structural  stresses  due  to  the  added  thickness,  but  an  additional  weight  penalfy  is 
imposed  upon  these  stronger  designs.  Because  DAPS3  does  not  consider  strength 
constraints,  the  smaller  L/OD  ratio  DAPS3  designed  shells  are  limited  by  stresses 
which  exceed  the  material  yieU  strength.  In  contrast,  the  THESIS  code  considers 
both  strength  and  buckling  constraints.  As  can  be  seen  in  Figures  14, 16,  18,  20,  22, 
and  24,  both  types  of  constraints  remain  satisfied  throughout  the  range  of  pressures 
for  all  L/OD  ratios.  The  THESIS  shell  weight  may  be  heavier,  but  the  design  is 
feasible  throughout  the  range  of  pressures. 
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B.  RING-STIFFENED  SHELL 


The  results  for  the  ring-stiffened  shell  are  presented  in  the  same  format  as  the 
monocoque  shell:  objective,  design  variables,  and  constraints.  There  are  now  three 
design  variables  (instead  of  one)  and  eight  constraints  (instead  of  six).  This  results 
in  a  more  complex  optimization  problem  since  the  design  space  is  now  three- 
dimensional  (due  to  three  design  variables).  The  simple  iterative  technique 
employed  by  DAPS3  produced  nearly  the  same  objective  (sheU  weight)  as  the 
THESIS  program,  which  used  a  numerical  optimization  technique.  However,  the 
optimum  design  variables  for  DAPS3  and  THESIS  differed.  Specifically,  the 
HEIGHT  and  WIDTH  of  the  ring-stifiener  varied  considerably.  A  graphical 
summary  of  the  results  are  seen  in  Figures  25  through  48. 

1.  Objective:  DAPS3  vs.  THESIS 

The  first  optimum  design  element  examined  is  the  objective:  optimum  shell 
weight.  Inspection  of  Figures  25  and  26  show  that  essentially  the  same  objective  was 
achieved  by  DAPS3  and  THESIS.  Some  slight  weight  reduction  is  seen  in  the  1:1 
and  3:1  THESIS  weight  curves,  which  can  be  attributed  to  a  lesser  weight  ring  design 
(that  is,  design  variables  of  HEIGHT  and  WIDTH).  This  will  be  discussed  in  the 
design  variable  results.  In  general,  it  can  be  concluded  that  the  iterative  design 
technique  essentially  achieves  the  same  objective  as  the  numerical  optimization 
technique.  The  other  design  elements  (design  variables  and  constraints)  must  be 
examined  to  reveal  if  numerical  optimization  provides  any  advantage. 


52 


WEIGHT  [LB/FT]  K  WEIGHT  [LB/FT] 


DAPS3  OPTIMUM  WEIGHT 


ONE  RING  SHELL  i/oo 


PSTR  [PSI] 


DAPS3  Optimum  Weight  (Rlng'StifTened  Shell) 


THESIS  OPTIMUM  WEIGHT 


ONE  RING  SHELL  uco 


PSTR  [PSI] 


Figure  26.  THESIS  Optimum  Weight  (Rlng'StifTened  Shell) 


2.  Design  Veriebles:  DAPS3  vs.  THESIS 

The  three  optimum  design  variables  produced  by  DAPS3  and  THESIS  are 
seen  in  Figures  27  through  32.  The  first  design  variable,  T,  is  compared  in  Figures 
27  and  28.  There  is  not  an  apparent  difference  between  the  two  programs.  Shell 
thickness,  T,  remains  a  monotonicalfy  increasing  variable  with  increasing  external 
hydrostatic  pressure.  Since  shell  thickness  is  the  primary  source  of  shell  weight,  the 
optimum  thickness  curves  look  very  similar  to  the  weight  curves,  as  was  seen  in  the 
monocoque  shell  case. 

The  addition  of  a  single  ring-stiffener  also  contributes  to  the  weight  curves, 
although  not  very  noticeably.  Within  the  surveyed  pressure  range  (500  psi  to  5800 
psi),  both  DAPS3  and  THESIS  design  a  relatively  small  single  ring-stiffener  into  the 
shell  design.  The  small  size  of  this  ring-stiffener,  in  relation  to  the  entire  shell,  is  the 
reason  that  the  weight  curves  resemble  only  the  shell  thickness  curves.  However,  this 
small  ring-stiffener  increases  the  allowable  external  hydrostatic  pressure  by  about  500 
psi.  That  is,  compared  to  a  similar  L/OD  ratio  monocoque  shell,  the  ring-stiffened 
shell  design  can  withstand  a  load  of  5800  psi  (instead  of  5300  psi  for  the  monocoque 
shell)  before  the  material  fails  by  yield. 

The  design  of  this  ring-stiffener  consists  of  optimizing  the  variables 
* ’EIGHT  and  WIDTH  such  that  the  shell  weight,  which  includes  the  weight  of  the 
ig,  is  minimized  and  all  constraints  are  satisfied.  As  seen  in  Figures  29  through  32, 
the  optimized  variables  of  ring  HEIGHT  and  ring  WIDTH  are  quite  different 
between  the  DAPS3  and  THESIS  codes. 
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T  [IN]  T  [IN] 


0  1000  2000  3000  4000  5000  6000 

PSTR  [PSI] 

DAPS3  Optimum  Shell  Hilckness  (T) 


THESIS  OPTIMUM  THICKNESS 


ONE  RING  SHELL  i/op 


ngure  28.  THESIS  Optimum  Shell  Thlcknesf  (T) 
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ngure  29.  DAPS3  Optimum  Ring  Height  (HEIGHT) 


I 


ngure  31.  DAPS3  Optiirium  Ring  Width  (WIDTH) 


Figure  32.  THESIS  Optimum  Ring  Width  (WIDTH) 
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Because  the  design  variables  HEIGHT  and  WIDTH  are  non*linear  in  the 
objective  function,  neither  variable  is  monotonically  increasing  or  decreasing  with 
respect  to  the  external  hydrostatic  pressure  load.  This  accounts  for  the  wide  anay 
of  ring  heights  and  widths  for  each  L/OD  ratio,  as  designed  by  the  THESIS  code 
(see  Figures  30  and  32).  By  contrast,  the  ring  heights  for  each  L/OD  ratio,  as 
designed  by  the  DAPS3  code,  are  generally  increasing  with  increasing  pressure  (see 
Figure  29),  and  the  ring  widths  for  each  L/OD  ratio  are  generally  decreasing  (see 
Figure  31).  Thus,  a  noticeable  difference  exists  in  the  way  that  the  two  programs 
optimize  non-linear  design  variables. 

At  this  point,  it  is  seen  that  the  same  objective  may  be  achieved  by 
different  design  variables.  The  selection  of  these  optimum  design  variables  depends 
on  the  particular  technique  employed:  iterative  (DAPS3)  or  numerical  optimization 
(THESIS).  As  in  the  case  of  monocoque  shell  design,  an  examination  of  constraint 
activity  by  the  two  different  programs  may  reveal  if  the  numerical  optimization 
technique  yields  any  advantages. 

3.  Constraints:  DAPS3  vs.  THESIS 

The  ring-stiffened  shell,  as  designed  by  THESIS,  has  eight  design 
constraints:  three  for  buckling,  four  for  strength,  and  one  for  geometry.  Recall  that 
DAPS3  designs  are  only  constrained  by  buckling  considerations.  This  difference 
results  in  noticeable  strength  constraint  violations  by  DAPS3  designs.  A  graphical 
comparison  of  all  eight  constraints  is  shown  in  Figures  33  through  48. 


DAPS3  G(l)  CONSTRAINT 

ONE  RING  SHELL 


Figure  33.  DAPS3  Constraint:  Symmetric  Local  Buckling  (PC) 


Figure  34.  THESIS  Constraint:  Symmetric  Local  Buckling  (PC) 
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Figure  35.  DAPS3  Constraiit.:  Non^Symmetric  Local  Buckling  (PIB) 


THESIS  G(2)  CONSTRAINT 


1000  2000  3000  4000  5000  6000 

PSTR  [PSl] 


THESIS  Constraint:  Non-Symmetrlc  Local  Buckling  (PIB) 


DAPS3  G(3)  CONSTRAINT 

ONE  RING  SHELL 


Figure  37.  DAPS3  Constraint:  Non^Symmetrlc  General  Buckling  (PGEN) 


Figure  38.  THESIS  Constraint*  Non'Symmetrlc  General  Buckling  (PGEN) 
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DAPS3  G(4)  CONSTRAINT 

ONE  RING  SHELL 


Figure  39.  DAPS3  Constraint:  Stress  (STR) 


Figure  40. 


THESIS  Constraint:  Stress  (STR) 


DAPS3  G(5)  CONSTRAINT 

ONE  RING  SHELL 


Figure  41.  DAPS3  Constraint:  Stress  (SIGVOF) 


THESIS  G(5)  CONSTRAINT 

ONE  RING  SHELL  i/oo 


Figure  42.  THESIS  Constraint:  Stress  (SIGVOF^ 
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DAPS3  G(7)  CONSTTRAINT 

ONE  RING  SHELL 


l/ao 


Figure  45. 


DAPS3  Constraint:  Stress  (SRF) 


THESIS  G(7)  CONSTRAINT 

ONE  RING  SHELL 


PSTR 


Figure  46.  THESIS  Constraint  Stress  (SRF) 
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DAPS3  G(8)  CONSTRAINT 

ONE  RING  SHELL 


Figure  47.  DAPS3  Constraint:  Ring  Aspect  Ratio 


THESIS  G(8)  CONSTRAINT 


PSTR 


Figure  48.  THESIS  Constraint:  Ring  Aspect  Ratio 


66 


An  inspection  of  Figures  33  through  38  shows  that  both  DAPS3  and 
THESIS  generally  satisfy  all  three  buckling  constraints.  A  visible  improvement  is 
seen  in  the  g,  constraint,  non-symmetiic  general  buckling.  A  comparison  of  Figure 
37  with  Figure  38  shows  that  this  constraint  is  better  adhered  to  by  the  THESIS 
code.  In  Figure  37,  the  DAPS3  g)  constraint  shows  slight  violations  (just  above  the 
zero  line)  but  is  generalfy  active  or  satisfied.  In  contrast.  Figure  38  shows  that  the 
THESIS  g3  constraint  is  always  active  or  satisfied.  Some  improvement  is  realized  by 
using  the  numerical  optimization  technique. 

Strength  constraints  are  considered  next.  These  are  shown  in  Figures  39 
through  46.  Notice  that  the  THESIS  g^(j  «  4, 5,6,7)  constraints  do  not  show  any 
violations:  all  designs  insure  that  stresses  do  not  exceed  the  material  yield  strength. 
In  contrast,  the  same  DAPS3  constraints  show  violations  in  two  cases.  The  first 
violation  is  seen  in  Figure  39,  which  is  the  g4  constraint  (stress  at  mid>bay).  This 
constraint  is  rather  small  in  comparison  with  the  g«  constraint  (stress  at  the  shell  and 
ring  junction  inner  fiber)  seen  in  Figure  43.  This  constraint  was  identified  as  the 
limiting  stress  location,  which  means  that  it  is  the  first  location  where  the  shell  would 
fail  by  yield.  At  higher  pressures,  the  g^  constraint  is  always  active  in  the  THESIS 
designs,  but  some  large  violations  occur  in  DAPS3  designs  at  these  same  pressures 
(see  Figure  43).  Thus,  use  of  a  numerical  optimization  technique,  which  easily 
incorporates  strength  constraints,  has  shown  an  advantage  over  the  DAPS3  technique 
which  does  not  concurrently  consider  strength  constraints  and  buckling  constraints. 
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By  coDsidering  both  strength  and  buckling  constraints,  the  THESIS  codes 
sometimes  designs  a  thicker  (and  heavier)  shell  than  the  DAPS3  code.  The  result 
is  that  the  strength  constraints  become  more  active  (since,  at  higher  pressures,  the 
material  yield  strength  is  the  limiting  design  factor)  and  some  buckling  constraints 
become  less  active.  An  example  of  this  is  seen  in  Figures  36  and  38,  which  show  the 
gj  and  g,  constraints  falling  below  the  zero  line  (beaming  less  active).  The 
correqK>nding  pressures  where  the  gj  and  g,  constraints  become  less  active  are  also 
the  pressures  where  all  the  strength  constraints,  g^(}  ^  4,S,6,7),  approach  the  zero 
line  (become  more  active).  This  is  seen  in  Figures  40,  42,  44,  and  46.  Thus,  the 
THESIS  code  may  design  a  slightly  heavier  shell,  but  this  shell  design  does  not 
violate  any  buckling  gi  strength  constraints. 

Geometry  is  the  last  type  of  constraint  considered.  The  gg  constraint  limits  the 
maximum  aspect  ratio  of  the  ring*stiffeoer  to  4:1.  This  was  done  to  prevent  design 
of  a  "flimsy"  ring.  This  constraint  is  shown  in  Figures  47  and  48.  Except  for  a  few 
of  the  smaller  L/OD  ratio  designs,  this  constraint  was  generally  not  active.  Initial 
ring'Stiffened  shell  designs  by  the  THESIS  code  showed  this  constraint  to  be  active 
at  higher  pressures,  however  these  earlier  designs  were  made  with  a  different 
material  (AISI  1060  steel,  instead  of  the  A1  alloy  6061'T6  used  in  this  study).  In 
summary,  this  constraint  was  not  active,  and  all  ring-stiffeners  could  be  considered 
"stable". 
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y.  CONCLUSIONS 


The  primaiy  goal  of  this  research  was  to  identify  any  improvements,  in 
monocoque  and  ling'Stiffeoed  shell  designs,  due  to  incoiporation  of  a  numerical 
optimization  technique  into  the  DAPS3  program.  The  design  elements  compared 
in  this  study  were:  1)  objective  (shell  weight),  2)  design  variables  (shell  thickness,  ring 
height,  and  ring  width),  3)  constraints  (buckling,  strength,  and  geometry). 

The  results  show  that  the  objective  (shell  weight)  was  not  further  minimized  by 
use  of  the  numerical  optimization  technique.  Both  DAPS3  and  THESIS  arrive  at 
about  the  same  optimum  shell  weight,  although  the  design  variable(s)  which 
contribute  to  this  shell  weight  may  differ. 

In  the  case  of  the  monocoque  shell,  it  was  found  that  the  design  variable,  shell 
thickness  (T),  was  about  the  same  for  both  DAPS3  and  THESIS.  This  was  due  to 
the  simple  one-variable  design  space  in  which  either  the  DAPS3  iterative  technique 
or  the  THESIS  numerical  optimization  technique  could  easily  achieve  the  single 
optimum  thickness  which  corresponds  to  the  minimum  monocoque  shell  weight.  In 
just  a  few  cases,  the  thickness  (T)  differed  between  the  two  programs.  In  these 
rases,  the  THESIS  design  variable  T  was  thicker  than  the  same  DAPS3  design 
variable  because  a  thicker  shell  was  necessary  to  satisfy  strength  constraints. 

In  the  case  of  the  ring-stiffened  shell,  the  design  variable  T  was  also  found  to 
be  about  the  same  for  both  DAPS3  and  THESIS  designs.  However,  the  ring- 
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stiffener  designed  by  the  two  programs  differed  in  ring  HEIGHT  and  ring  WIDTH. 
This  difference  in  design  variables  is  not  signiEcant  with  respect  to  minimizing  the 
shell  weight.  That  is,  the  optimum  design  variables  produced  by  DAPS3  and 
THESIS,  although  different,  still  result  in  the  same  shell  weights.  In  the  more 
complex  three-variable  design  ^ace  of  this  q>timization  problem,  it  was  found  that 
HEIGHT  and  WIDTH  are  nonlinear  variables  in  the  objective  function.  This 
accounted  for  the  "up  and  down"  range  of  values,  seen  in  Figures  30  and  32,  for 
HEIGHT  and  WIDTH.  Although  HEIGHT  and  WIDTH  do  not  significantly 
contribute  to  shell  weight  (due  to  the  relatively  small  size  of  the  ring-stiffener),  these 
design  variables  made  a  difference  in  complying  with  the  design  constraints. 

The  primaiy  advantage  of  using  a  numerical  optimization  technique  for  shell 
design  was  realized  by  an  examination  of  the  constraints.  In  general,  both  DAFS3 
and  THESIS  satisfied  the  buckling  constraints  for  both  monocoque  and  ring-stiffened 
shells.  However,  a  significant  difference  was  seen  in  the  strength  constraints. 

In  the  monocoque  shell  case,  it  was  found  that  some  DAPS3  strength 
constraints  (especially  at  the  ring  and  shell  junction  inner  fiber  location)  are  violated 
(see  Figure  21).  In  contrast,  the  same  THESIS  strength  constraints  were  satisfied 
due  to  an  increase  in  shell  thickness,  as  discussed  above. 

In  the  ring-stiffened  shell  case,  it  was  found  that  similar  strength  constraint 
violations  occurred  in  DAPS3  designs.  These  constraint  violations  were  absent  from 
THESIS  designs.  Since  THESIS  shell  designs  consider  strength  constraints,  while 
DAPS3  designs  do  not,  the  THESIS  code  was  able  to  optimize  the  design  variables 
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(especially  ring  HEIGHT  and  ring  WIDTH)  such  that  there  were  no  strength 
constraint  violations. 

In  summary,  incorporation  of  a  numerical  optimnation  technique  into  DAPS3 
proved  to  be  a  worthwhile  venture.  The  primary  advantage  realized  is  that 
constraints  can  easily  be  added  and  tailored  to  fit  a  more  constrained  shell  design 
problem.  This  is  simply  done  by  adding  a  mathematical  constraint  function,  gj,  to  the 
optimization  model.  DAPS3  is  an  excellent  shell  design  program,  but  it  scnnewhat 
lacks  user  flexibility  in  specifying  additional  constraints  (such  as  strength  constraints). 
Another  advantage  realized  is  that  additional  design  variables,  such  as  shell  length 
(L)  and  shell  outside  diameter  (OD),  may  be  optimized  when  the  flexible  numerical 
optimization  technique  is  used. 

The  THESIS  code  was  not  written  to  refine  the  solutions  produced  by  DAPS3, 
but  rather  to  supplement  the  number  of  solutions  available.  In  this  age  of  high¬ 
speed  digital  computer  technology,  the  ease  of  employing  numerical  optimization 
techniques  could  be  the  key  to  discovering  alternative,  possibly  more  effective, 
solutions  to  traditional  problems. 
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VI.  RECOMMENDATIONS 


The  THESIS  code  was  a  first  endeavor  into  incorporating  numerical 
optimization  with  the  DAPS3  program.  Several  possibilities  exist  for  more 
sophisticated  optimization  research.  The  following  are  suggested  for  any  successive 
research: 


1.  New  design  variables.  Investigate  optimum:  number  of  rings,  shell  diameter, 
or  shell  length. 

2.  New  constraints.  One  example  is  maximum  allowable  shell  deflection. 

3.  New  materials.  Composite  (orthotropic)  materials  and/or  hybrid  construction 
materials.  Dramatic  shell  weight  reductions  may  be  realized  by  use  of  these 
materials.  The  DAPS3  "analysis”  and  "design”  modes  fully  supports  research 
of  these  materials. 

4.  New  optimization  str:  ^egies.  Several  optimization  strategies  exist  in  the  ADS 
program  which  may  improve  upon  the  SLP  strategy  used  in  this  study. 
Additionally,  a  parametric  optimization  technique  such  as  the  homotopy 
method  using  an  envelope  function  may  be  investigated  [Ref.  7]. 
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APPENDIX  A:  THESIS  CODE  LISTING 


This  appendix  provides  a  complete  listing  of  the  THESIS  code.  It  is  written 
in  FORTRAN  to  be  compatible  with  the  DAPS3  code,  which  is  also  written  in 
FORTRAN.  Also  included  is  the  modified  DAPS3  code,  which  is  used  as  a 
subroutine  (DAPS3  SUB)  by  THESIS. 
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C  . . . . . 

C  THESIS 

C  (AN  INTEGRATED  FORTRAN  CODE  FOR  OPTIMUM  DESIGN  OF  ISOTROPIC 
C  MONOCOQUE  AND  RING-STIFFENED  CIRCULAR  CYUNDRICAL  SHELLS 

C  SUBJECT  TO  EXTERNAL  HYDROSTATIC  PRESSURE  USING  A  NUMERICAL 

C  OPTIMIZATION  TECHNIQUE) 

. . . . . . . 

C 

C  THIS  FORTRAN  CODE  INTEGRATES  THE  PROGRAMS: 

C  I)  DESIGN  AND  ANALYSIS  OF  PLASTIC  SHELLS,  VERSION  3  (DAPS3) 

C  BYJ.R.RENZI 

C  2)  AUTOMATED  DESIGN  SYNTHESIS,  VERSION  1.10  (ADS) 

C  BY  G.N.  VANDERPLAATS 


C  PROGRAM  (1)  IS  A  SHELL  ANALYSIS/DESIGN  PROGRAM  WHICH  IS  USFD  TO 
C  CALCULATE  THE  BUCKLING  AND  STRENGTH  TERMS  CONTAINED  IN  THE 
C  DESIGN  CONSTRAINTS,  G(I). 

C  PROGRAM  (2)  RECEIVES  THE  TERMS  EVALUATED  BY  PROGRAM  (1)  AND 
C  EVALUATES  THE  DESIGN  CONST  RAINTS  AS  PART  OF  THE  OVERALL 
C  NUMERICAL  OPTIMIZATION  PROCESS. 

C  THE  OBJECnVE  OF  THIS  NUl  .lERICAL  OPTIMIZATION  PROCESS  IS  TO 
C  DESIGN  AN  OPTIMUM  (MI>!iMUM  WEIGHT)  RING-STIFFENED  OR 
C  MONOCOQUE  SHELL. 


C  - 

PROGRAM  THESIS 

INTEGER  NPL0T,NCURVE,NDESIGN.0RTH0,IPRINT,ISTRAT,10PT,I0NED 
REAL  ITHICK,IHE1GHT,IWIDTH 

OPEN  (20,FILE=’/DAPS3  !NP  Al’.STATUS-’OLD') 

READ  (20,’(4(15))’)  NPLOT,NCURVE,NDESIGN,ORTHO 
IF  ((NPLOT.EQ.O).AND.(NDESIGN.EQ.-1)) 

:  READ  (20,’(//////F10.3)’)  AR 
IF  ((NPLOT.EQ.O).AND.(NDESIGN.NE.-1)) 

;  READ  (20,’(/////F10.3)’)  AR 
IF  ((NPLOT.NE.O).AND.(NDESIGN.EQ.-1)) 

:  READ  (20,’(/////F10.3)’)  AR 
IF  ((NPLOT.NE.O).AND.(NDESIGN.NE.-1)) 

:  READ  (20,’(////F10J)’)  AR 

CLOSE(20,STATUS= ’KEtr’) 

C  INPUT  ADS  OPTIMIZATION  PARAMETERS: 

IPRINT  «  2020 
ISTRATe  6 
lOPT  -5 
lONED  -  6 


C  INPUT  INITIAL  DESIGN  FOR  OPTIMIZATION  SEARCH: 

IF  (AR.LE.0)  THEN 
ITHICK  -  1.1589 

CALL  MONO  (IPRINT.ISTRAT,IOPT,IONED,nHICIC) 

ELSE 

ITHICK  ■  0.191 
IHEIOHT  -  0.191 
IWIDTH  -  1.730 

CALL  RING  (IPRINT.ISTRAT,IOPT,IONED,ITHICK,IHEIGHT,IWIDTH) 
ENDIF 

END 


C - 

SUBROUTINE  MONO  (IPRINT,ISrRAT,IOPT,IONED,ITHICK) 
C - 


REAL  UOD.T4’C;PIB,R.RO,RHO.RHOF^YIUP^OSl,MOS2,MOS 
REAL  MSS1.MSS2,MSS3.MSS4.MSS,LS.1R, ITHICK 

INTEGER  NPLOT,NCURVE,r4DESIGN,ORTHO,NSAVE,ISEC,NXPNT,NYPNT,NFPNT 
INTEGER  IPRINT.ISTRAT,IOPT,IONED 
CHARACTER'70  TITLE 
PARAMETER  (PI  •  3.141593) 

DIMENSION  X(2l).  VLB(21),  VrjB(21),  G(IOO),  IDG(IOO),  IC(30). 

DF(I00),  A^2130),  WK(IOOOO),  IWK{5000) 

COMMON  AR,IR.YR,Lr.LS,ARM,l,RO,T,GXY3X^YEF,VXY,VYX,VPLAS,R, 

:  PSTR.STR,SIGX,SIGY,SIOXUIGYl^EIGHT,WIDTH,OD 

COMMON  /COUNT/  nERAT,INDIC,INDIC2,INDlC3, 

:  INDIC5,NDESIGN, ORTHO, NOKG.NCURVE,OK,nERAM,OK;4,nERAG,OKO,ISEC, 

:  NPLOT,NOKM,NRINGS.NjaLAM,INDIC7 


COMMON  EMOD/  STRESX(10),STRINX(10),STRESY(IO),STRINY(10). 

:  STRESF(!0),STRiNF(10),RATIOXRATIOY,NXPNT,NYPNT,NFPNT 

COMMON /GRAPK/ET(28),SIG(28),PB(28),PS(28),ES(28)A^A(6), 

•.  3BB(6),CCC(8),DDD(6)EEE{5)4^;T(7).0GG(4)J1HH(7), 000(8), ER(28). 
:  SIGF(28),SHAPEO(5),SHAPEI(3),WX(10I) 

COMMON  /LAYER/  A11,A22D11,D22 

COMMON  /STRESS/  STRX,STRY,PIB,NSAVEI,PGEN,NSAVEO 

COMMON  /HENRY/  SRF,WM,WF 

OPEN  (8.FILE- ’/THESIS  OUT  Al’,STATUS-'OLD’) 
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C  COMMENCE  READING  INPUT  VARIABLES  FROM  (204nLE-'/DAPS3INPAr): 
OPEN  (204^ILE-'/DAPSJ  INP  Al’.STATUS-’OlD’) 

READ  (20400)  NPLOT.NCURVE4IDESION.OR1HO 
READ  (2040)  TTILE 

WIUTE(d,’(A/A)’)  ‘COMMENCING  OPTIMIZATION  OF:’,  TITLE 
IF  (NPLOT.EQX)) 

:  READ  (2040)  HHH 
IF  (NDESION.EQ.-1) 

:  READ  (2041)  A11,A224>11.D22 
READ  (2041)  EX,VYX,VPLAS,SI0X,0XY3Y£F 
READ  (2043)  PSTR,SIOY,SIGX1.SIOY1  J>EUSTRX.STRY,NSAVE 
READ  (2041)  L,OD,T.RHO,RHOF 
IF  (PELLEO.) 

:  READ  (2041)  AR.1R,YR.LF.LS,ARM 
IF  (NDESIGN.GT.0) 

:  READ  (2043)  RINGS.SECDEP.DWJ)ASP,RINGSI4lINOSF4>RINOS,ISEC 
IF  (NCURVE.EQ.4)  GOTO  2 

READ  (2042)  SIGXNOM.SlGYNOM.NXPNT,NYPNT,NFPNT 
READ  (2045)  (STRESX(I),  I  -  1,NXPNT) 

READ  (2045)  (STRINX(I),  I  -  l.NXPNT) 

IF  (NCURVE.EQ.1)  GOTO  2 

READ  (2045)  (STRESY(I).  I  -  1 WNT) 

READ  (2045)  (STRINY(I),  1  «  1,NYPNT) 

IF  (NCURVEfQ.2)  GOTO  2 

READ  (2045)  (STRESF(I),  I  -  l.NFPNT) 

READ  (2045)  (STRINF(I),  I  -  144FPNT) 

500  FORMAT  (BZ,10I5) 

50  FORMAT  (B2,A) 

51  FORMAT  (BZ,8F104) 

52  FORMAT  (BZ,2F10J4I5) 

53  FORMAT  (BZ,7F104,I5) 

55  FORMAT  (BZ,10F8D) 

2  CONTINUE 

CIjOSE(20,STATUS-’KEEP’) 

NRA  -  21 
NCOLA-  30 
NRWK  -  lOOOO 
NRIWK  -  5000 

IGRAD  -  0 
NDV  •  1 
NCON  -  6 

C  THE  LOWER  BOUND  IS  THE  MIN  ALLOWED  SHELL  THICKNESS  (INCHES) 
VLB(1)»  0.10 
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C  THE  UPPER  BOUND  IS  SET  TO  1^18  IN.  so  THAT  THE  MAX  ALLOWED 
C  T/R  RATIO  IS  20%.  (THIS  VALUE  IS  VALID  FOR  THE  FIXED  VALUE 
C  10  IN.  OUTSIDE  SHELL  DIAMETER  US'  3  IN  THIS  STUDY.) 

VUB(l) "  1.818 

DO  3  I  «  1.6 
IDO(I)-0 
3  CONTINUE 

WRITE  (8,'(A)')  • - ’ 

WRITE  (8,’(5XAn.T20AIl.T32An)*) 

:  TSTRAT  -  *4STRAT,TOPT  -  ’.lOPT.TONED  -  MONED 
WRITE  (8,’(A)')  ’ - ’ 

X(l)-  riHICK 

WRITE  (8,’(/1XAT20AF7J.T50A/)’)  'INITIALIZED  AT;’, 

:  ’THICKNESS(IN)  »  ’.XCl), ’(DESIGN  VARIABLE  #1)’ 

WRITE  (8.’(lXA'nOAF7.2,2XAF7.2/)’) 

:  BOUNDARY  CONSTRAINTS—’. 

:  ’VARIABLE  #1:’,  VLB(I).  TX)’,  VUB(l) 

INFO  -  -2 

CALL  ADS  (INFO, ISTRAT.IOPT,IONED,IPRINT,lGRAD,NDV,NCON,X, VLB, 

vub,obj,g,idg,ngt,iqdfanra,ncoiawk,nrwk,iwk,nriwk) 

IWK(2)  -  0 

10  CALL  ADS  (INF0,ISTRAT,10PT,I0NED.IPRINT,IGRAD,NDV.NC0N,X.VLB. 

VUB,0BJ,G.IDG,NGT,IC43FANRA,NC0IAWK,NRWK,IWK,NRIWK) 
if  (INFO  .EQ.  0)  GOTO  20 

T  -  X(l) 

RO  •  OD/2. 

R  ■  RO  -  T/2. 

CALL  DAPS3  (PC,RHO,RHOF,nTLE,SIGVOF,SIGVIF,USEVOL,CL) 

OBJ  -  (24.*PI/L)  •  (R*T»L*BHO  +  (RO-YR)*AR*(L/LF4.)*RHOF) 

C  THE  0(1)  CONSTRAINT  INSURES  AXISYMMETRIC  COLLAPSE  PRESSURE  IS 
C  LESS  THAN  ANALYSIS  PRESSURE 
G(l)  -  1.  -  (PC/PSTR) 

C  THE  0(2)  CONSTRAINT  INSURES  ELACTIC  MONOCOQUE  INSTABUJl Y 
C  OR  INTER-BAY  INSTABILITY  PRESSURE  IS  LESS  THAN  ANALYSIS  PRESSURE 
0(2)  ■  1.  -  (PIB/PSTR) 

C  THE  0(3)  CONSTRAINT  INSURES  MID-BAY  MEMBRANE  STRESS  LEVEL 
C  (VON  MISES  STRESS,  FROM  SIOXM  AND  S101M,  PER  YIELD)  IS  LESS 
C  THAN  100%  (OF  MATERIAL  YIELD  SIHENOTH  SIOX) 

G(3)  "  (STR/100.)  - 1. 
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C  THE  0(4)  CONSTRAINT  INSURES  VON  MISES  S1RESS  IN  MEMBER  AT 
C  OUTSIDE  FRAME  IS  LESS  THAN  100%  OF  MATERIAL  YIELD  STRENGTH 
0(4)  -  (SIGVOF/100.)  - 1. 

C  THE  G(S)  CONSTRAINT  INSURES  VON  MISES  STRESS  IN  MEMBER  AT 
C  INSIDE  FRAME  IS  LESS  THAN  100%  OF  MATERIAL  YIELD  STRENGTH 
G(5)  -  (SIGVnVlOO.)  - 1. 

C  THE  0(6)  CONSTRAINT  INSURES  THAT  HOOP  STRESS  IN  THE  FRAME  (RING) 
C  IS  LESS  THAN  MATERIAL  YIELD  STRENGTH 
IF  (ORTHO^jO)  then 
0(6)  -  ABS(SRF/SIGX)  •  1. 

ELSEIF  ((OR'IHOmi)AND.(SRFXTjO))  THEN 
0(6)  -  ABS(SRF/SIGY)  - 1. 

ELSEIF  ((ORTHO.EQ.l)AND.(SRF.aTi)))  THEN 
G(6)  -  (SRF/SIGYI)  - 1. 

ENDIF 
GOTO  10 
20  CONTINUE 

C  DEFINE  LIMITING  MARGIN  OF  STABILTrY  (MOS) 

MOSl  -  PC/PSTR 
MOS2  -  PIB/PSTR 
MOS  «  MIN(MOSl.MOS2) 

C  DEFINE  LIMITING  MARGIN  OF  STRENGTH  (MSS) 

MSSl  -  STR/100. 

MSS2  -  SIOVOF/100. 

MSS3  «  SIOVIF/IOO. 

IF  (ORTHO.EQ^)  THEN 
MSS4  «  ABS(SRF/SIGX) 

ELSEIF  ((ORTHO.EQ.1)AND.(SRF.LT.O))THEN 
MSS4  -  ABS(SRF/SIGY) 

ELSEIF  ((ORTHO.EQ.l)J^ND.(SRF.GT^))THEN 
MSS4  -  SRP/SIGYI 
ENDIF 

MSS  -  MAX(MSSl.MSS2^SS3>fSS4) 

WRITE  (8,’(1X^H***  OPTIMIZATION  OF:/IX^)*)  TITLE 
WRITE  (8,’(T6,A,F7A1XA0’)  ’®  ANALYSIS  PRESSURE  -  '4»SrR,’PSr 
IF  (ORTHO.BQi))  THEN 

WRITE  (S/Cre^)*)  'ISOTROPIC  MATERIAL  PROPERTIES;' 

ELSE 

WRITE  (8,'(T6,A)’)  'ORTHOTROPIC  MATERIAL  PROPERTIES:’ 

ENDIF 

WRITE  (8,’(T11^)')  ’»>SHELL’ 

IF  (ORTHO.EQ.O)  THEN 
WRITE  (8,'(T14A2PE9^IX^)’) 

:  'EX  »  ’EX,TSr 
ELSE 

WRITE  (8,'(2(TI4A2PE92,1X^,T14.AE9X1X^)’) 
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:  'EX  -  ’^TSr, 

:  'EY  -  ’3Y,TSr, 

:  TEF  -  'JEF.TSr 

ENDIF 

WRITE  (8,’(T14A1PE92,1X^)’) 

:  ’OXY  -  ’.GXY.’PSr 
IF  (ORlHO.EQjO)  THEN 
WRITE  (8,’(T14AF8.0,1X^)’) 

:  'YffilD  -  ’,SIOX,TSr 
ELSE 

WRITE  (8,’(T14A3(T14APW.lXA),T14AF8.aiXA)') 

:  ’YIELD:’, 

:  ’AXIAL(C)-  ’,SIGX.’PSr, 

:  ’AXIALCO-  ’,SIGXl,TSr, 

:  ’HOOP  (Q-  ’.SIGY.TSr, 

:  ’HOOP  (I)-  ’.SIGYLTST 
ENDIF 

WRITE  (8,’(T14AFJJ.1X,A)’) 

:  ’DENSITY  -  ’.RHO.’LB/CU.IN.’ 

WRITE  (8.’(T14AF6.1,1X,A)’) 

:  ’  -  ’,RHO*(12**3),’LB/CU.FT.’ 

WRITE  (8.’(3(T14AF63/))*) 

:  ’VXY  -  ’.VXY, 

:  ’VYX  -  ’.VYX. 

:  'VPLAS  -  ’  VPLAS 

WRITE  (8,’(3(ixA,F94/).lXAP9.4,T50A2(lX,A4J9^/))’) 

:  ’CRITICAL  LENGTHCIN)  -  CL, 

:  ’LENGTH(IN)  -  L, 

:  ’OUT.DIAMETER(IN)  -  OD, 

:  ’THICKNESS(IN)  -  X(l), ’(OPTIMIZED  VARIABLE  #1)’, 

:  ’RING  SPACING(IN)  »  LF, 

;  ’INTERBAY  WIDTHCIN)  -  LS 
WRITE  (8,’(2<1XAF6A,IXA))’) 

:  ’SHELL  RADIAL  DEFLECHON  @  MIDBAY  -  ’,  WM,  ’IN’, 

©FRAME  - ’,WF, ’IN’ 

IF  (MOS.EQ.MOS1)  THEN 
WRITE  (8,’(1XAF6.4)’) 

’LIMITING  MARGIN  OF  STABILITY  (PC/PSTR)  -  ’,  MOS 
ELSEIF  (MOS.EQ.MOS2)  THEN 
WRITE  (8.’(1XAF6.4)’) 

’UMITTNO  MARGIN  OF  STABILITY  (FIByPSTR)  -  ’,  MOS 
ENDIF 

IF  (MSS.EQ.MSSI)  THEN 
WRITE  (8,’(1XAW-4/)’) 

’UMITTNO  MARGIN  OF  STRENGTH  (STR/YIELD)  -  ’,  MSS 
ELSEIF  (MSS.EQ.MSS2)  THEN 
WRITE  (8,’(1XAF6.4/)’) 

’UMITTNO  MARGIN  OF  STRENGTH  (SIOVOF/YIELD)  -  MSS 
ELSEIF  (MSS.EQ.MSS3)  THEN 
WRITE  (8,’(1XAF6.4/)’) 
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’UMITINO  MARGIN  OP  CTRENOTH  (SIOVIP/YIELD)  *  MSS 
ELSEIF  ((MSS.EQ.MSS4).AND.(ORTHO£Qi)))  THEN 
WRITE  (8,’(1XAF6.4/)’) 

UMITING  MARGIN  OF  STRENGTH  (SRF/SIGX)  >  \  MSS 
ELSEIF  ((MSS.EQ.MSS4)j\ND.(ORTHO£Q.l).AND.(SRFXT.O))  THEN 
WRITE  (8,’(1XAF6.4/)’) 

’UMITING  MARGIN  OF  STRENGTH  (SRF/SIGY)  -  MSS 
ELSEIF  ((MSS£QA1SS4)  AND.(ORlHO£Q.l)  JKND.(SRF.GT.O))  THEN 
WRITE  (8,’(lXAF6.4/)‘) 

’UMITINO  MARGIN  OF  STRENGTH  (SRF/SIGYl)  -  MSS 
ENDIF 

WRrre  (8,’(6(1XA1PE10J,T25A/))*) 

:  ’G(l)  »  ’.G(1),’AXISYMMETRIC  COLLAPSE  PRESSURE  :  PC, 

:  ’G(2)  •  ’.G(2),'MONOOOQUE/IN1ERBAY  INSTABIUrY  PRESURE  :  PIB’, 

:  ’G(3)  -  ‘,G(3),’VON  MISES  MIDBAY  MEMBRANE  mESS  :  STR’. 

:  ’0(4)  «  ’,0(4), 'VON  MISES  OUTSIDE  FRAME  MEMBRANE  STRESS;  SIGVOF, 

:  ’G(S)  -  ’,G(5),'VON  MISES  INSIDE  FRAME  MEMBRANE  STRESS  :  SIOVIF, 

:  ’0(6)  -  ’.G(6),’HOOP  STRESS  IN  THE  FRAME  :  SRF 

WRITE(8,’(2(1XA.F10J/),1XAF103,T50A/1XAF10J//1XAI6//)’) 

:  ’INTERNAL  VOLUME(CU.FTyFT.)  -  ’,  USEVOL, 

:  ’INTERNAL  VOLUME(CU.FT)  -  ’,  USEVOL*L/12, 

:  ’WEIGHT(LB/FT)  «  ',  OBJ,’(OPTTMIZED  OBJECTIVE)’, 

:  ’WEIGHT(LB)  -  ’,  OBJ*L/12, 

:  ’ADS  ITERATIONS  -  ’,  IWK(28) 

IF  (LLT.CL)  THEN 
WRITE  (8,’(3(1X,A/))’) 

;  ’NOTE:  SHELL  LENGTH  IS  LESS  THAN  MINIMUM  CRITICAL  LENGTH.’, 

:  ’END  PLATE  EFFECTS  MAY  AFFECT  REPRESENTATIVE  STATES  OF  STRESS’, 
:  ’ATMED-SPANOFTHESHELL.’ 

ENDIF 

CLOSE  (8,STATUS-'KEEF) 

END 


C - 

SUBROUTINE  RING  (IPRINT,ISTRAT,IOPT,IONED,nHICK,IHEIGHT,IWIDTH) 
C - 


REAL  L,0D,TJ»C,PIB,R,R0,RH0,RH0FAR,YILLF>I0S1,M0S2,M0S3,M0S 

REAL  MSS1>ISS2,MSS3,MSS434SS,LS,1R.ITHICK,IHEIGHT,IWIDTH 

INTEGER  NFLOT,NCURVEJ4DESION,ORTHO,NSAVE,ISEC,NXPNT,NYPNT,NFPNT 

INTEGER  N0R,IPRINT,ISTRAT,I0PT,I0NED 

character*70  title 

PARAMETER  (PI  -  3.141593) 

DIMENSION  X(21),  VLB(21),  VUB(21).  0(100),  IDO(IOO),  IC(30), 
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DF(IOO),  A(21^),  'VK(10000),  1WK(5000) 


COMMON  ar,ir,yr,lf,ls,arm,l,ro,t,gxy^ey^f.vxy,vyx,vplas,r, 

:PSTR,CTR.SlGX,SIGY,SIGXl,SIGYl^lGHT,WIDTH,OD 
COMMON  /COUNT/  ITCRAT.INDIC.INDIC2,INDIC3, 

:INDIC5,NDESIGN.0RTH0J40K0.NCURVE,0ICITERAM,0KM.ITERAG,0KG.ISEC; 
:  NPLOT,NOKM^RINOS,NEGL\M.INDIC7 

COMMON  /EMOD/  CTRESX(10),STRINX(10),CTRESY(10),STRINY(10), 

:  STRESF(10),STRINF(10),RAT!OX.RATIOY4aXPNT,NYPNT,NFPNT 

COMMON /GRAPH/ET(28).SIG(28),PB(28),PS(28),ES(28)AAA(6), 
:BBB(6),OCC(8),DDD(6)^EE(5)JTF(7).GGG(4)3HH(7), 000(g), ER(28). 

:  SIGP(2^,SHAPEG(S}.SHAPEI(3),WX(101) 

COMMON  /LAYER/  A11,A224)11,D22 

COMMON  /STRESS/  STRX.STRY,PIB,NSAVEI,PGEN,NSAVEG 

COMMON  /HENRY/  SRF.WM.WF 

OPEN  (84?ILE- '/THESIS  OUT  Al*. STATUS- ’OLD  ) 

C  COMMENCE  READING  INPUT  VARIABLES  FROM  (204TLE-’/DAPS3  INP  Al*): 
OPEN  (204^ILE-’/DAPS3  INP  Al'.STATUS-’OLD’) 

READ  (20400)  NPLOT,NCURVENDESION,ORTHO 
READ  (2040)  TITLE 

WRITE(d,’(A/A)’)  'COMMENCING  OPTIMIZATION  OF:’,  TITLE 
IF  (NPLOT.EQ.O) 

:  READ  (2040)  HHH 
IF  (NDESIGN.EQ.-1) 

:  READ  (2041)  A11,A22,D11.D22 
READ  (2041)  EX,VYX,VPLAS,SIGX,GXY,EYEF 
READ  (2043)  PSTR,SIGY,SIGXl,SIGYl,PEL,STRX,SrRY,NSAVE 
READ  (2041)  L,OD,T,RHO,RHOF 
IF  (PELLE.0.) 

:  READ  (2041)  AR,IR,YR,LF,LS, ARM 
IF  (NDESIGN.GT.0) 

:  READ  (2043)  RINGS, SECDEP.DW4)ASP,RINGS1,RINGSFJ>RINGS,ISEC 
IF  (NCURVE.EQ.4)  GOTO  2 

READ  (2042)  SIGXNOM,SlGYNOM,NXPNT,NYPNT,NFPNT 
READ  (2045)  (STRESX(I),  I  -  1,NXPNT) 

READ  (2045)  (STRINX(I),  I  -  1,NXPNT) 

IF  (NCURVE.EQ.1)  GOTO  2 

READ  (2045)  (STRESY(I),  I  -  14IYPNT) 

RPAD  (2045)  (STRINY(I),  I  -  1,NYPNT) 

IF  (NCURVE.EQ.2)  GOTO  2 

READ  (2045)  (STRESF(I).  I  -  1,NFPNT; 
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READ  (20^5)  (CTRINF(I).  I  -  1,NFPNT) 

500  FORMAT  (BZ,10IS) 

50  FORMAT  (BZ»A) 

51  FORMAT  (BZ.8F103) 

52  FORMAT  (B42F1033I5) 

53  FORMAT  (BZ,7F103,I5) 

55  FORMAT  (B410F8.0) 

2  CONTINUE 

CLOSE(20,STATUS-'KEEF) 

NRA  -21 
NCOLA-30 
NRWK  -  10000 
NRIWK  -  5000 

IGRAD  -  0 
NDV  -3 
NCON  -  8 

C  THE  LOWER  BOUNDS  ON  (1)  THICKNESS,  (2)  RING  HEIGHT,  AND 
C  (3)  RING  WIDTH  ARE  ESTABLISHED  (INCHES) 

VLB(1)-  0.10 
VLB(2)-  0.05 
VLB(3)-  0.125 

C  THE  UPPER  BOUNDS  ON  (1)  THICKNESS,  (2)  RING  HEIGHT,  AND 
C  (3)  RING  WIDTH  ARE  ESTABUSHED  (INCHES) 

VUB(l)  -  1.818 
VUB(2)-  10. 

VUB(3)-  5. 

DO  3  I  -  1,8 
IDG(I)  -  0 

3  CONTINUE 

WRITE  (8,’(A)’)  ’ - ’ 

WRITE  (8,’(5XAI1.T20AI1.T32AI1)’) 

:  TSTRAT  -  MSrRAT,TOPT  -  MOPT,TONED  -  ’,10NED 
WRITE  (8,’(A)’)  ’ - ' 

X(l)-  nHICK 
X(2)-  WEIGHT 
X(3)-  rWIDTH 

WRITE  (8,’(/lXA3(T20AF73.T50A))’)  ’INITIALIZED  AT:’, 

:  ’THICKNESS(IN)  -  ’,X(1), ’(DESIGN  VARIABLE  #1)’, 

:  ’RING  HEIOHT(IN)  -  ’,X(2),’(DESION  VARIABLE  #2)’, 
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:  'RING  WIDlll(IN)  - ’.X(3V(DES10N  VARIABLE  #3)’ 
WRITE  (8,'(lX^(T30AF72,2XAF7Jl/)y) 

:  'BOUNDARY CONSTRAINTS—’, 

;  'VARIABLE  #1:’,  VLB(1),  TO',  VUB(l), 

;  'VARIABLE  #2:’,  VLB(2),  TO',  VUB(2), 

:  'VARIABLE  #3:’.  VLB(3),  TO’,  VUB(3) 


INFO  -  -2 

CALL  ADS  (INF04STRAT40PT40NED4PRINT,I0RAD,NDV,NC0N,X,VLB, 
VUB.OBI,04D044aT4C;DFANRA4ICOLA,WK4WWK,IWK,NRIWK) 
IWK(2)-0 

10  CALL  ADS  (INF04STRAT40PT40NED4PRINT,I0RAD,NDV,NC0N.X.VLB. 

VUB,0BJ,0,ID044aT,IC4>FANRA4<C0LA,WK44RWK,IWK,NRIWK) 

IF  (INFO  .EQ.  0)  GOTO  20 

T  -  X(l) 

HEIGHT  -  X(2) 

WIDTH  -  X(3) 

RO  -  ODy2. 

R  -  RO  -  T/2. 

YR  -  T  +  HEIGHTy2. 

AR  -  HEIGHT  •  WIDTH 

CALL  DAPS3  (PCRHO.RHOF,TnLE,SIGVOF,SIGVlF,USEVOL,CL) 

OBJ  -  (24 •  (R'TL'RHO  +  (RO.YR)*AR*(lVLF-l.)*RHOF) 

C  THE  G(l)  CONSTRAINT  INSURES  AXISYMMETRIC  COLLAPSE  PRESSURE  IS 
C  LESS  THAN  ANALYSIS  PRESSURE 
G(l)  -  1. .  (PQPSTR) 

C  THE  0(2)  CONSTRAINT  INSURES  ELASTIC  MONOCOQUE  INSTABILITY 
C  OR  INTER-BAY  INSTABILmr  PRESSURE  IS  LESS  THAN  ANALYSIS  PRESSURE 
G(2)  -  1.  -  (PIB/PSTR) 

C  THE  G(3)  CONSTRAINT  INSURES  THAT  EGI  PRESSURE  (RING-STIFPENED 
C  SHELL)  IS  LESS  THAN  ANALYSIS  PRESSURE 
0(3)  -  1.  -  (PGEN/PSTR) 

C  THE  G(4)  CONSTRAINT  INSURES  MID-BAY  MEMBRANE  STRESS  LEVEL 
C  (VON  MISES  STRESS,  FROM  SIOXM  AND  SIOTM,  PER  SIGX)  IS  LESS 
C  THAN  100%  (OF  MATERIAL  YIELD  STRENGTH  SIGX) 

0(4)  -  (SfR/100.)  - 1. 

C  THE  0(3)  CONSTRAINT  INSURES  VON  MISES  STRESS  IN  MEMBER  AT 
C  OUTSIDE  FRAME  IS  LESS  THAN  100%  OF  MATERIAL  YIELD  STRENGTH 
0(5)  -  (SIOVOF/100.)  -  1. 

C  THE  0(6)  CONSTRAINT  INSURES  VON  MISES  STRESS  IN  MEMBER  AT 
C  INSIDE  FRAME  IS  LESS  THAN  100%  OF  MATERIAL  YIELD  STRENGTH 
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G(6)  -  (SIOVIF/100.)  -  1. 


C  THE  G(7)  CONSTRAINT  INSURES  THAT  HOOP  STRESS  IN  THE  FRAME  (RING) 
C  IS  LESS  THAN  MAIBRIALYIELO  STRENGTH 
IF  (ORTHO.EQ.O)  THEN 
G(7)  ■  ABS(SRF/SIOX)  -  1. 

ELSEIF  ((OR'IHO.EQ.l)j\ND.(SRF.LT.O))THEN 
0(7)  «  ABS(SRF/SIOY)  •  1. 

ELSEIF  ((OR'IHO.EQ.l)AND.(SRF.GT^))THEN 
G(7)  -  (SRF/SIGYl)  -  1. 

ENDIF 

C  THE  0(8)  CONSTRAINT  LIMITS  RING  ASPECT  TO  4:1  HEIGHT/WIDTH  RATIO 
G(8)  -  (HEIGHT/WIDTH)/4..  1. 

GOTO  10 
20  CONTINUE 

C  DEFINE  LIMITING  MARGIN  OF  STABILITy  (MOS) 

MOSl  »  PC/PSTR 
MOS2  -  PIB/PSTR 
MOS3  -  PGEN/PSTR 
MOS  B  MIN(MOSl.MOS2,MOS3) 

C  DEFINE  LIMITING  MARGIN  OF  STRENGTH  (MSS) 

MSSl  -  STR/100. 

MSS2  •  SIGVOF/100. 

MSS3  «  SIGVIF/100. 

IF  (OR1HO.EQ.O)  THEN 
MSS4  >  ABS(SRF/SIGX) 

ELSEIF  ((OR’mO.EQ.l)AND.(SRF.LT.O))THEN 
MSS4  -  ABS(SRF/SIGY) 

ELSEIF  ((ORTHO,EQ.1).AND.(SRF.GT.O))THEN 
MSS4  -  SRF/SlGYl 
ENDIF 

MSS  >  MAX(MSS1,MSS2.MSS3,MSS4) 

C  DEFINE  NUMBER  OF  RINGS  (NOR) 

NOR  -  INT((L/LF)-1.) 

WRITE  (8,'(1X^H***  OPTIMIZATION  OF:/lX,A)’)  TTIIE 
WRITE  (8,’(T6yvF7.0,lXAO’) '®  ANALYSIS  PRESSURE  -  ’4»STR,TSr 
IF  (ORTHO.EQ.O)  THEN 

WRITE  (8,'(T6,A)’)  'ISOTROPIC  MATERIAL  PROPERTIES:’ 

ELSE 

WRITE  (8,’(T6,A)’)  ’ORTHOTROPIC  MATERIAL  PROPERTIES:’ 

ENDIF 

WRITE  (8,’(T1 1  AT43^)’)  ’>  >  >SHELL:’,’>  >  >  RING:’ 

WRITE  (8,’(T14A2PE9J,1XA,T46,A^9^1X,A)’) 

:  ’EX  -  ’JEX’PSr,  'EF  -  ’EF.’PSI' 

IF  (ORTHO.EQ.1) 

: WRITE  (8,’Cn4A2PE9^1X,A)’)  ’EY  -  ’^Y.TSl’ 
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WRITE  (8,’Cn4AlPE9^1XAT46AE9^1X.A)’) 

:  ’OXY  •  '.GXY.TSr,  'OXY  -  ’.OXY/PSI’ 

IF  (ORTOO.EQ.O)  THEN 

WRITE  (8,’(TI4AF8.0,1XAT46AF8.0,1XA)’) 

’YIELD  «  ’.SlOX'PSr, ’YIELD  ■  ’^IOX,’PSr 
ELSE 

WRITE  (8,’(T14AT46A)’) 

•YIELD:*, •YIELD:’ 

WRITE  (8,’(T14AF8.0,1XAT46AF8.0.IXA 
T14AF8X),1XAT46AFSj0.1XA)') 

’HOOP(C)  -’,SIOY,*PSr.*HOOP(0  -’.SIOY.TSr, 

•HOOPCI)  -*,SIOYl,*PSI’,’HOOPCn  -’.SIOYl.’PSr 
WRITE  (8,’(T14AF8.0,lXAn4AF8.0,lXA)’) 
•AXIAL(q»’,SIOX.’PSr,’AXIAL(T)-’.SIOXl.*PSI’ 

ENDBP 

WRITE  (8,’(T14AF5  J.1XAT46AFS  J.1XA)’) 

:  ’DENSITY-  ’.RHO.TB/CU.IN.’, 

:  ’DENSITY-  ’.RHOF.’LB/CU.IN.’ 

WRITE  (8,’(T14AF6.1,1XAT46AF6.1,1XA)*) 

:  ’  -  ’.RHO*(12**3).’LB/CUJT.’, 

;  ’  -  ’,RHOF*(12**3),’LB/CU.Fr.’ 

WRITE  (8,’(3(T14AF6J,T46AF6.3/))’) 

:  ’VXY  -  ’.VXY,  'VXY  -  ’.VXY, 

:  ’VYX  -  ’.VYX,  ’VYX  -  ’.VYX, 

:  ’VPIAS  -  ’.VPLAS,  ’VPLAS  -  ’.VPLAS 
WRITE  (8,’(3(1XAF9.4/)3(1X,A.F9.4,T50A), 
1XAF9.4,1X,1H(,I2,1X6HRINGS)/1XAF9.4/)’) 

:  ’CRITICAL  LENG'IH(IN)  -  CL, 

:  ’LENGTH(IN)  -  L, 

:  'OUT.DIAMETER(IN)  -  OD. 

:  ’THICKNESS(IN)  -  ’.XCl), ’(OPTIMIZED  VARIABLE  #1)’, 

:  ’RING  HEIGHT(IN)  -  ’,X(2),’(OPTIMIZED  VARIABLE  #2)’, 

:  ’RING  WIDTH(IN)  -  ’,X(3), ’(OPTIMIZED  VARIABLE  #3)’, 

:  ’RING  SPACING(IN)  «  ’aJ,NOR, 

:  ’INTERBAY  WIDTH(IN)  -  ’,LF  •  X(3) 

C  NOTE:  INTERBAY  WIDTH  HAS  BEEN  REDEFINED  AS  LS  -  LF  •  X(3) 

C  SO  THAT  THE  ORIGINAL  DAPS3  INP  VALUE  OF  LS  IS  NOT  USED. 
WRITE  (8,’(2(1XAF6.4,1XA))’) 

:  ’SHELL  RADIAL  DEFLECTION  @  MIDBAY  -  ’,  WM,  ’IN’, 

:  ’  @  FRAME  -  ’,  WF,  ’IN’ 

IF  (MOS.EQ.MOS1)  THEN 
WRITE  (8,’(1XAF6.4)’) 

’LIMITING  MARGIN  OF  STABILITY  (PC/PSTO)  -  ’,  MOS 
ELSEIF  (MOS.EQ.MOS2)  THEN 
WRITE  (8,’(1XAF6.4)’) 

’UMITTNO  MARGIN  OF  STABILnY  (PIB/PSTR)  -  ’,  MOS 
ELSEIF  (M0S.EQ.M0S3)  THEN 
WRITE  (8,’(1XAF6.4)') 

’UMITING  MARGIN  OF  STABILnY  (POEN/PSTR)  -  ’,  MOS 
ENDIF 
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IF  (i.^o6.EQ.MSSl)  THEN 
WRITE  (8,’(1XAF6^’) 

'UMITINO  MARGIN  OF  SIHENGTH  (STR/YIELD)  -  MSS 
ELSEIF  (MSS.EQJ4SS2)  THEN 
WRITE  (8,’(1XAF6.4/)*) 

’UMITING  MARGIN  OF  STRENGTH  (SIGVOF/YIELD)  -  \  MSS 
ELSEIF  (MSS.EQ.MSS3)  THEN 
WRITE  (8,*(1XAF6.4/)’) 

’UMmNG  MARGIN  OF  CTRENOTH  (SIOVIF/YIELD)  -  MSS 
ELSEIF  ((MSS.EQ.MSS4).AND.(ORTHO£QX)))TTIEN 
WRITE  (8,’(1XAF«.4/)') 

‘UMITINO  MARGIN  OF  STRENGTH  (SRF/YIELD)  •  MSS 
ELSEIF  ((MSS.EQ.MSS4).AND.(ORTHOEQ.1).AND.(SRFXT.O))  THEN 
WRITE  (8.'(1XAF8.4/)’) 

TJMmNO  MARGIN  OF  ^IRENGTH  (SRF/SIGY)  «  MSS 
ELSEIF  ((MSS.EQ.MSS4)  AND.(ORTHO£Q.1).AND.(SRF.OT.O))  THEN 
WRITE  (8,’(1XAF6.4/)*) 

•UMITING  MARGIN  OF  STRENGTH  (SRF/SIGYl)  -  MSS 
ENDIF 

WRITE  (8,’(8(1XA1PE10J,T25A/))’) 

:  ’G(l)  »  ’.G(1).*AXISYMMETRIC  COLLAPSE  PRESSURE  :  PC, 

:  ’G(2)  -  •,G(2),'MONOCOQUE/INTERBAY  INSTABILITY  PRESURE  :  PIB*, 

:  ’G(3)  ■  ’.G(3),'ELASnC  GENERAL  INSTABILITY  PRESSURE  :  PGEN’, 

:  *0(4)  «  •.0(4).*V0N  MISES  MIDBAY  MEMBRANE  STRESS  :  STR', 

:  ’G(5)  »  ’.G(5).’VON  MISES  OUTSIDE  FRAME  MEMBRANE  STRESS:  SIGVOF, 

:  *0(6)  «  *,0(6).* VON  MISES  INSIDE  FRAME  MEMBRANE  STRESS  :  SIGVIF, 

:  ’G(7)  «  ’,G(7).’HOOP  STRESS  IN  THE  FRAME  :  SRF, 

:  ’0(8)  -  ',0(8), 'RING  ASPECT  RATIO(HT/WDT)  .LE.  4:1’ 
WRITE(8,’(2(1XAF10J/),1XAF10,3,T50,A/1XAF103//1X,A,I$//)’) 

:  'INTERNAL  V0LUME(CU.Fr4T.)  -  USEVOL 
:  'INTERNAL  VOLUME(CU.FT)  -  USEVOL*U12„ 

:  ’WEIGHT(LB/FT)  -  OBJ,'(OPTIMIZED  OBJECITVE)’, 

:  ’WEIGHT(LB)  -  OBJ*U12, 

:  'ADS  ITERATIONS  -  IWK(28) 

IF  (LLT.CL)  THEN 
WRITE  (8,’(3(1XA))’) 

:  NOTE:  SHELL  LENGTH  IS  LESS  THAN  MINIMUM  CRITICAL  LENGTH.’, 

:  ’END  PLATE  EFFECTS  MAY  AFFECI  REPRESENTATIVE  STATES  OF  STRESS’, 
:  AT  MID  SPAN  OF  THE  SHELL’ 

ENDIF 

CLOSE  (8,STATUS-’KEEF) 

END 
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C  . . . . . . . 

SUBROUTINE  DAPS3  (PC,RH0,RH0F,Tni^S10V0F,SIGVIF,USEV0L,CX) 

C  . . . . . . . 

C  (THIS  SUBROUTINE  IS  ESSENT7  \LLY  THE  SAME  AS  THE  ORIGINAL  DAPS3 
C  PROGRAM.  SOME  MODIFICATIONS  WERE  NECESSARY  TO  INTEGRATE  THIS 
C  CODE  WITH  THE  THESIS  CODE.  SUCH  MODIFICATIONS  ARE  NOTED  IN 
C  COMMENT  LINES.) 

DIMENSION  RSAVE(99),WSAVE(99) 

INTEGER  OK,OKM.OKG,ORTHO.RCOUNT,OKR 
REAL  IR,ULS,LF 
CHARACTER‘70  TnLE(I) 

COMMON  AR,IR.YR,LF,LSARM.L,RO.T.GXYEX,EYEF,VXY,VYX,VPLAS,R, 
:PSTR,STR,SIGX,SIGY,SIGXl,SIGY14IEIOHT,WIDTH,OD 
COMMON  /COUNT/  nERAT,INDICINDIC2,INDIC3, 

:  INDIC5,NDESIGN, ORTHO, NOKG.NCURVE,OIC,lTERAM,OKM,nERAG.OKG,ISEC, 
:  NPLOT, NORM, NRINGS.NEGLAM, INDICT 

COMMON  /EMOD/  STRESX(10),STRINX(10),STRESY(10),STRINY(10), 
:STRESF(10).STRINF(10),RATIOX,RATIOY,NXPNT,NYPNT,NFPNT 
COMMON  /GRAPH/ET(28),SIG(28).PB(28),PS(28),ES(28)AAA(6), 

:  BBB(6),0CC(8),DDD(6),EEE(5),FFF(7).GGG(4),HHH(7),000(8),ER(28), 

:  SIGF(28).SHAPEG(5).SHAPEI(3).WX(101) 

COMMON  /LAYER/  A1 1^224)11.022 

COMMON  /STRESS/  STRX,STRY,PIB,NSAVE1J*GEN,NSAVEG 

COMMON  /HENRY/  SRF.WM.WF 

PI  «  4.*ATAN(1.) 

ITERAT  «=  0 
RINGSI  *  0. 

RINGSF  »  0. 

DRINGS  -  0. 

3  IF  (RINGSI.GT.O.)  GO  TO  4 

C  UNIT®  6,  FILE  ®/THESIS  D/>iPSOUT  Tl:  IS  DEFINED  BY  THE  THESIS  EXEC 
C  UPON  EXECUTION  OF  THE  MAIN  PROGRAM  (THESIS  FORTRAN  Al) 

100  CONTINUE 
10000  CONTINUE 

4  IF  (NDESIGN.EQ2)  GO  TO  5 
WRITE  (6,60)  nTLE(l) 

IF  (NDESIGN.EQ.l)  WRITE  (6,64) 

5  IF  (NDESIGN.EQ.2)  WRITE  (6,60)  'nTLE(l) 

IF  ^DESIGN.EQ.2)  WRITE  (5,645) 

OK  -  -1 
ITERAM-  0 
ITERAG-  0 
OKM®  0 
OKG-  0 
RCOUNT- 0 
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NOKM- 0 
NOKG  s  0 
NEOLAM  «  0 
PEL-0. 

IF  (RINGSI.GT.O.)  GO  TO  6 

NXPNT-  0 

NYPNT  -  0 

NFPNT  »  0 

DO  7  I  -  1,10 

STRESX(I)-  0. 

SrrRINX(I)-  0. 

STRESY(I)-  0. 

STRINY^-  0. 

STRESF(I)»  0. 

STRINF(I)»  0. 

7  COmiNUE 

C  AR  IS  COMMENTED  OUT  TO  PREVENT  CHANGING  INPUT  VALUE  FM  DAPS3  INP 
C  AR  =  0. 

RINGS  «  0. 

10007  IF  (NCURVE.EQ.4)  GO  TO  6 

10014  IF  (NCURVE.GT.1  .AND.  NYPNT.EQ.O)  WRITE  (6,6100) 

IF  (NCURVE.EQ.3  .AND.  NFPNT.EO.0)  WRITE  (6,6110) 

6  CONTINUE 

IF  (PEL.GT.0..AND.NDES1GN.NE.0)  WRITE  (6,6000) 

IF  (PEL.GT.0..AND.NDESIGN.NE.0)  GO  TO  999 
IF  (NDESIGN.EQ.2)  WRITE  (6,646)  PSTR 

IF  (NCURVE.EQ.2)  EF  -  EY 

IF  (ORTHO-EQ.!)  GO  TO  8 

EY  -  EX 

EF  -  EX 

SlGYl  -  SIGX 

SIGXl  -  SIGX 

SIGY  -  SIGX 

STRY  -  STRX 

SIGYNOM  -  SIGXNOM 

GXY  -  EX/(2.*(1.+VYX)) 

IF  (VPLAS.GT.  0.5)  VPLAS  -  0.5 

8  CONTINUE 

IF  (ORTHO.EQ.O  AND.  NCURVE.EQ.3)  EF  -  EFSAVE 
IF  (ORTHO.EQ.O  AND.  NOIRVE  LT^)  RHOF  -  RHO 
STRY  -  ABS(STRY) 

STRX  -  ABS(SniX) 

C  IF  PEL  .OT.  0.,  THEN  AT  THIS  POINT  BOTH  STRX  AND  STRY  MUST  BE 
C  GREATER  THAN  ZERO. 

IF  (EX.LT.10.)  EX  -  10. 
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IF  (EY.LT.10.)  EY  -  10. 

VXY  •  VYX*EX/EY 

IF  (NDESION.EQ.-1)  VXY  -  VYX*A22/A11 
IF  (NDESIGN.EQ..1)  GO  TO  14 
IF  (VXY.GT.SQRT(EX/EY)/2.)  VXY  -  SQRT(EX/EY)/2. 

VYX  -  VXY*EY/EX 

IF  (VYX.OT.SQRT(EY/EX)/2.)  VYX  -  SQRT(EY/EX)/2. 

VXY  •  VYX*EX/EY 
14  CONTINUE 

IF  (NDESIGN.GT.-l)  GO  TO  16 

IF  (VXY.GT.SQRT(A22/All)/2.)  VXY  -  SQRT(A22/A11)/1 
VYX  -  VXY*A11/A22 

IF  (VYX.OT.SQRT(All/A22)/2.)  VYX  «  SQRT(A11/A22)/1 
VXY  •  VYX*A22/A11 

16  CONTINUE 

IF  (VPLASiE.O.)  VPLAS  •  VYX 
IF  (SIGXNOM.LE.O.)  SIGXNOM  -  SIGX 
IF  (SIGYNOM.LE.O.)  SIGYNOM  »  SIGY 
RATIOX  -  SIGX/SIGXNOM 
RATIOY  -  SIGY/SIGYNOM 
DO  9  I  -  1.40 
WSAVE(I)  «  0. 

9  RSAVE(I)  -  0. 

IF  (NDESIGN.GT.0  .AND.  PSTR.LE.O.)  WRITE  (6,6300) 

IF  (NDESIGN.GT.0  .AND.  PSTR.LE.0.)  GO  TO  997 
RINGS  -  AINT(RINGS) 

IF  (NDESIGN.LT.2)  GO  TO  10 

C  REMEMBER,  IN  NDESIGN  -  2  OPTION.  NDESIGN  IS  RESET  FROM  2  TO  1 
C  AFTER  SEARCHING.  RMAX  IS  THE  MAXIMUM  NUMBER  OF  RINGS  CONSIDERED. 
IRING  -  0 

IF  (RINGSI.GT.0.  .AND.RINGSF.OT.0.  .AND.(RINGSF-RINOSI).OT.O.) 

;  IRING  -  1 
RINGS  -  RINGSI 
RMAX  -  RINGSF 
IF  (IRING.EQ.1)  GO  TO  17 
RINGS  -  1. 

RINGSI  -  1. 

RMAX  «  AINT(L/1.5) 

IF  (RMAXGT.66.)  RMAX  -  66. 

RINGSF  -  RMAX 

17  CONTINUE 

DRINGS  -  AINT(RMAX/I5.) 

IF  (DRINOS.LT.3,)  DRINGS  -  3. 

IF  ((RINGSF-RINOSI).LT,3.  .AND.IRING.EQ.1)  DRINGS  -  1. 

10  CONTINUE 

C  IF  (NDESION.GT.0)  PRINT  *,  *  CALL  RIBBIT  FROM  DAPS3  NEAR  LABEL  10* 

IF  (NDESIGN.OT.0)  CALL  RIBBIT(RINOS,SECDEPJ)W,DASP) 

1  OK  -  OK  -«•  1 

IF  (NCURVE.EQ.4)  OK  -  1 
INDIC  ■  0 
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INDIO  «  0 
INDIO  -  0 
INDICT  -  0 
RO  -  ODy2. 

R  -  RO  -  m 
RR  RO  •  YR 

IF  (NDESION.EQ..1)  EXOEY  -  A22/A11 
IF  (NDESIGN.GT.-l)  EXOEY  -  EX/EY 
CL  •  PI*(R**2  •  T**2  /(3.*(1.  -  VXY*VYX))*EXOEY)**25 
IF  (AR.LE.0.)  THEN 
LS-  L 
LF-L 
END  IF 

IF  (LF.GT.0.  .AND.  NDESIGN.LE.O)  RINGS  -  L/LF  -  1. 

IF  (RINGS.GT.0.  .AND.  (RINGS-AINT(RINGS)).GT.J)  RINGS  »  RINGS  +1. 

RINGS  -  AINT(RINGS) 

NRINGS  «  RINGS 
IF  (PffrR.LE.O.)  OK  -  1 

IF  (NDESIGN.GT.0  .AND.  OKEQ.O)  INDIO  -  1 
IF  (INDIO.EQ.l)  GO  TO  15 
INDICT  -  1 

C  SETTING  INDICT  -  1  ALLOWS  THE  FULL  FIITH  ORDER  SET  OF  EQUATIONS  TO 
C  BE  USED  IN  GENERL  IN  THE  FINAL  SET  OF  CALCULATIONS  AFTER  THE 
C  DESIGN  ITERATIONS. 

IF  (NDESIGN.NE..1)  WRITE  (6,61)  EXEY,EF,GXY,VXY,VYX,VPLAS 
IF  (NDESION.EQ.-1)  WRITE  (6,612)  A11,A224>11,D22,GXY,EF,VXY,VYX 
IF  (PELGT.0.)  GO  TO  11 
VOLF  -  0, 

IF  (YR.GT.0.)  VOLF  -  2.*PI*(RO-YR)*AR*(ULF.l.) 

VOLS  -  2.*PI*R*T*L 

WEIGHT  -  24.*P1/L  •  (R'TL'RHO  +  (RO-YR)*AR*(LLF-l.)*RHOF) 

USEVOL  -  (L*PI*OD**2/4.  •  (VOLF+VOLS))/(LM44.) 

WRITE  (6,615)  WEIGHT,USEVOL,RHO,RHOF 

1 1  IF  (PELGT.0.)  WRITE  (6,605) 

IF  (PELGTX).)  GO  TO  12 

IF  (AR,OT.O.)  WRITE  (6,62)  ODL,CLTLF 
IF  (AR.LE.0,)  WRITE  (6,67)  ODL,CLT 
IF  (AR.GTX).)  WRITE  (6,625)  LS 
IF  (AR,OTi).)  WRITE  (6,63)  IR,AR,RR 

12  CONTINUE 

IF  (RINOS.OT.O..AND.NDESION.OTX))  WRITE  (6,65)  HEIOHT,WIDTH,NRINGS 

SECnON  -  T  +  HEIGHT 

IF  (WIDTH.GT,0.)  ASPECT  -  HEIGHT/WIDTH 

IF  (RINOS.GT.O,  .AND.  NDESION.OT.O)  WRITE  (6,66)  SECnON^VSPECT 

IF  (PELLE.0.)  GO  TO  13 

WRITE  (6,672)  ODX,CLTJ»EL 

WRITE  (6,673)  SrRY,STRX,PSTR 

GO  TO  40 

13  CALL  RENSA1JP(PC,0,SIOVOF,SIGVIF) 

IF  (PSTR.LE.0.)  GO  TO  997 
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15  CONTINUE 
INDIC  -  1 

C  PROGRAM  CONTROL  (GENERAL,  INTERBAY,  AND  AXISYMMETRIC  BUCKLING.) 
IF(AR.LE.O.)  GO  70  70 
INDIO  -  I 

C  PRINT  *,  •  CALL  RENKEN  FROM  DAPS3  NEAR  LABEL  70" 

CALL  RENK£N(PINELG) 

IF  (NCURVE.E0.4)  WRITE  (6,655)  PGEN.NSAVEG 
IF  (NCURVE  .EO.  4)  WRITE  (6,662)  SHAPEG 
IF  (NCURVE.EQ.4  .AND.  NEGLAM.EQ.I)  WRITE  (6,663) 

70  CONTINUE 

C  PRINT  *,  ’  CALL  RENVON  FROM  DAPS3  NEAR  LABEL  70  " 

CALL  RENVON(PINELB) 

IF  (NCURVE.EQ.4.ANDAR.GT.0.)  WRITE  (6,692)  PIB,NSAVEI 
IF  (NCURVE.EQ.4.ANDAR.LE.0.)  WRITE  (6,693)  PIBJ4SAVEI 
IF  (NCURVE  .EQ.  4)  WRITE  (6,694)  SHAPEI 
C  IF  (NCURVE.E0.4)  GO  TO  100 
40  CONTINUE 

IF  (INDIO.EQ.O)  OKG  «  1 
1  -  0 

IF  (INDICEQ.0)  CALL  REDUCE(PEL,NSAVE,I) 

IF  (IND1C.EO.O)  GO  TO  997 
IF  (INDIC5.NE.1)  GO  TO  90 

C  FROM  HERE  TO  90  IS  ONLY  SEEN  IN  A  DESIGN  OPTION. 

OKR  «  1 

CALL  RENSALP(PC,l,SlGVOF.SIGVIF) 

PDUM  -  AMINl(PCPINELB) 

CALL  rrERM(PDUM) 

IF  (OKM.EQ.l  .AND.  INDIO.EQ.O)  GO  TO  90 
IF  (INDIO.EQ.l)  CALL  ITERG(PINELG,SECDEP.DW,DASP) 

OKM  -  0 
OKG  -  0 

IF  (INDIO.EQ.O)  OKG  -  1 

75  ITERAT  -  ITERAT  +  1 

CALL  RENSALP(PC,l,SIGVIF,SIGVOF) 

C  PRINT  •,  "  CALL  RENVON  FROM  DAPS3  NEAR  LABEL  75  " 

CALL  RENVON(PINELB) 

PDUM  -  AMINI(PC,PINELB) 

CALL  ITERM(PDUM) 

IF  (NOKM.E0.2)  WRITE  (6,6200) 

IF  (NOKM.EQ.2)  THEN 
OKM  -  I 
OKG  «  1 
END  IF 

IF  (OKM.EQ.O)  GO  TO  75 
IF  (INDIO.EQ.O)  GO  TO  90 

76  ITERAT  -  ITERAT  +  I 

C  PRINT  ■  CALL  RENKEN  FROM  DAPS3  NEAR  LABEL  76" 

CALL  RENKEN(PINELG) 

CALL  ITERG(PINELG,SECDEP,DW.DASP) 
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IF  (NOKG.EQ.  2)  WRITE  (6,6200) 

IF  (NOKG.EQ.  2)  THEN 
OKM  B  1 
OKG-  1 
END  IF 

IF  (OKG.EO.O)  THEN 
OKM  -  0 
ITERAM  -  0 
END  IF 

IF  (OKG.EQ.O)  GO  TO  76 
IF  (OKM.EQ.O)  THEN 
OKG  «  0 
ITERAG  «  0 
END  IF 

IF  (OKM.EQ.O)  GO  TO  75 

C  WHEN  NDESIGN  -  1  (WHICH  INCLUDES  THE  FINAL  DESIGN  OF  NDESIGN  »  2 
C  WHERE  NDESIGN  WAS  RESET  TO  1).  INDIC7  SERVES  THE  SOLE  PURPOSE  OF 
C  CAUSING  FIRST  THE  SHELL  TO  BE  DESIGNED  USING  CUBIC  AND  THEN 
C  REFINED  FURTHER  TO  THE  FINAL  OPTIMUM  DIMENSIONS  USING  EIGZF  (SEE 
C  GENERL).  THEREFORE  ALL  FINAL  DESIGNS  USE  EIGZF. 

INDIC7  ■  INDIC7  +  1 

IF  (INDIC7.EQ.1  .AND.  NOKG.LT.2  .AND.  NOKM.LT.2)  THEN 
OKM  e  0 
OKG  «  0 
ITERAM  «  0 
ITERAG  =  0 
END  IF 

IF  (OKM.EQ.O  .AND.  NDESIGN.EQ.l)  GO  TO  75 
C  INDIC7  DOES  NOTHING  WHILE  NDESIGN  IS  STILL  2  . 

C  PRINT  DAPS3  ",RO,T,AR.YR,HEIGHT.WIDTH,LS 
IF  (NDESIGN.EQ.l)  GO  TO  90 

C  SEARCH  FOR  OPTIMUM  WEIGHT  (NDESIGN  «  2)  IS  CONTROLLED  HERE. 
OKM-  0 
ITERAM-  0 
OKG-  0 
ITERAG-  0 
NOKM- 0 
NOKG  -  0 

RCOUNT  -  RCOUNT  +  I 
IF  (RCOUNT.EQ.99)  GO  TO  997 

WEIGHT  -  24.*PI/L  •  (R*T*L*RHO  +  (R0-YR)*AR*(L;LF1.)*RH01) 
RSAVE(RCOUNT)  -  RINGS 
WSAVE(RCOUNT)  -  WEIGHT 
C  A«SECOND()  -  CALCULATTS  CPU  SECONDS 
WRITE  (6,666)  RINGS,WEIGHT,ITERAT,A 
84  RINGS  -  RINGS  DRINGS 
NRINGS  «  RINGS 
IF  (RINGS.GT.RMAX)  GO  TO  85 
CALL  RIBBIT(RINGS,SECDEP,DW,DASP) 

GO  TO  75 
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85  CXJNTINUE 

CALLRCLOSE(RSAVE,WSAVE,RINOS,WElGHT,RINGSI,DRINGS,RMAX, 

1  RCOUNT.OKR) 

IF  (OKR.EQ.l)  GO  TO  84 
ITERAT-  0 
NDESIGN  »  1 
RINGSI  -  RINGS 
RINGSF  «  RINGS 
DRINGS  -  1. 

WRINGS  *=  RINGS 
WRITE  (6,6400)  NRINGS 
G0T03 
90  CX)NTINUE 
NOKM  *  0 
NOKG  -  0 

C  FOR  AN  ANALYSIS  OF  A  MONOCOQUE  SHELL,  OKG  -  1  HERE  SINCE  INDIO  = 
C  (ABOVE  AFTER  40)  AND  OK  «  0  (FIRST  TIME  THRU)  BUT  OKM  NE.  1  SO 
C  THE  PROGRAM  STILL  DOES  NOT  GO  TO  I  . 

IF  (OKM.EQ.l.AND.OKG.EO.l.AND.OK.EQ.0)GO  TO  1 

60  FORMAT  (//5H . ,A70) 

605  FORMAT  (IHO) 

61  FORMAT  ( /5X.‘El ASOC  MODULI  V  7X.‘EX  -  *,2PE1033X,*EY  «  *, 

:  ElOJ,  3X,’EF  «  '.ElOJ,  3X.’GXY  =  MPE103  /  5X,‘POISSON  RATIO 
:S*  ,3X,*VXY  «  ’,  0PF53,  3X,*VYX  «  '4?5J3X.‘VPLAS  -  ’J^5  J) 

612  FORMAT  (  /5X‘A11  -  MPEI0J,3X’A22  -  ’,E10J3X,’Dn  -  ’,E103, 

:  3X,’D22  «  ’,E103/ 5X*GXY  «  MPEI0J,3X.’EF  -  ',E1033X, 

:  ’VXY  *  ’.0PF5J3X,‘VYX  »  ’,F5 J) 

615  FORMAT  (5X.'SHELL  WEIGHT  •  ’,F72.’  LB/FT .  INTERNAL  VOLUME  » 

:  F8.4,’  CU.FTVFT.75X,’(SKIN  DENSITY  -  ',F6J,7X, 

:  ’RING  DENSITY  -  *  ,F6.3.*)7 ) 

62  FORMAT(5X,’OUTSIDE  DIAMETER  -  LENGTH  ■  *J^83,*  CRITI 

:CAL  LENGTH  «  •,F53y5X,*BAY  THICKNESS  -  ',F53,’  RING  SPACING  - 
:’,F6.3) 

625  FORMAT  (5X,’UNSUPPORTED  BAY  WIDTH  -  ’T^63,2X,*(USED  FOR  INTERBAY 
:BUCKLING)’/) 

63  FORMAT  (5X,’RING  PROPERTIES’  /8X,’IX  =  •,1PE1033X,  ’AREA  -  ’, 

:  E1033X,  ’RADIUS  TO  CG.  -  ’,0PF73) 

64  FORMAT(y21X,’THIS  SHELL  IS  DESIGNED  BY  THE  PROGRAM.  ’) 

645  FORMAT  (/4X,’THE  DAPS3  SEARCH  OPTION  FOR  THE  OPTIMUM  NUMBER  OF  RIN 
:GS  GOES  AS  FOLLOWS.  ) 

646  FORMAT  (  25X,’(DESIGN  PRESSURE  -  ’,F6.0,’  PSI)’///  8X, 

.  ’NUMBER  OF  RINGS’, 5X, 'WEIGHT  (LB/FT)’, 5X.’ITERAT10NS’3X, 

:’CP  SECONDS’/) 

65  FORMAT  (8X.’RING  HEIGHT  «  ’,F73,’  RING  WIDTH  -  ’,F73,4X,’(’. 

:  13,’  RINGS)’) 

655  FORMAT  ( /  5X,’ELASTIC  GENERAL  INSTABILITY  PRESSURE  -  ’,IPE103, 
:3X,’MODE  -  ’,12) 

66  FORMAT(8X,’SECnON  DEPTH  -  ’,F73,’  ASPECT  RATIO  -  ’,F73) 

662  FORMAT  (  5X,’A1«’,1PEI03,’  BI»’,EI03,’  CI-’,EI03,’  B2«’, 

:  E103,’  C2»’,E103) 
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663  FORMAT  (  5X,’(A  NEGATIVE  EIGENVALUE  WAS  IGNORED.)’  ) 

666  FORMAT  (I3X,F4.0.13X,F83,I0X.I7,6X,F10.3) 

67  FORMAT(5X,’OUTSlDE  DIAMETER  »  LENGTH  »  ’,F83,’  CRITI 

:CAL  LENGTH  -  ’,F53/5X,*SKIN  THICKNESS  »  •,F53/5X, ’(ANALYSIS  IS  F 
;OR  A  MONOCOQUE  SHELL.)’) 

672  FORMAT(5X,  OUTSIDE  DIAMETER  -  ',F63,*  LENGTH  -  ’4^83,’  CRITI 
:CAL  LENGTH  =  ’,F53/5X,’SKIN  THICKNESS  -  •,F53//SX,’THE  ELASTIC  B 
:UCKLING  PRESSURE  (PRECALCULATED)  *  •.IPE103.’  PSI.  ) 

673  FORMAT  (/5X,’THE  PRECALCULATED  STRESS  LEVELS  IN  THE  SHELL  «  ’.IPEl 
:0.3,’  PSI  (HOOP)  AND  ’/SX.1PE10.3.’  PSI  (AXIAL)  FOR  AN  APPUED  PR 
:ESSURE  OF  ’3103,’  PSI.  ’/) 

692  FORMAT(  /SX’ELASTIC  INTER-BAY  INSTABILITY  PRESSURE  «  ’.1PE103, 
:3X,’MODE  =  ’.12) 

693  FORMAT  ( y5X.’ELASrnC  MONOCOQUE  SHELL  INSTABILITY  PRESSURE  «  ’, 

;  1PE103,2X,’MODE  -  ’,12) 

694  FORMAT  (5X,’BUCKLED  SHAPE  ...  ’,2X‘A  »  ’,1PE103,2X,’  B  « 

:  E103,2X.’  C  «  ’,0PF53) 

699  FORMAT  (/5X.’(THIS  SHELL  REQUIRED  ’,17,’  NUMERICAL  ITERATIONS.)’) 

6000  FORMAT  (/3X,’USER  WARNING  -  YOU  CANNOT  USE  THE  DESIGN  OPTION  WHEN 
:  PEL  IS  GREATER  THAN  0.  ’/) 

6100  FORMAT  (  5X,’NOTE  -  NO  POINTS  ON  THE  Y-CURVE  WERE  READ  IN  SINCE  N 
:YPNT  »  0  .’) 

6110  FORMAT  (  5X,’NOTE  -  NO  POINTS  ON  THE  F-CURVE  WERE  READ  IN  SINCE  N 
:FPNT  »  0  .’) 

6200  FORMAT  (  5X.’*  ATTEMPTS  AT  CONVERGENCE  FOR  THIS  DESIGN  ARE  DISCON 
:TINUED.  ’) 

6300  FORMAT  (/SX.  USER  MESSAGE  -  NDESIGN  MUST  BE  0  OR  -1  WHEN  PSTR  IS 
:NEGAT1VE.’) 

6400  FORMAT  ( //  8X,’ ALTHOUGH  ’.13,’  RINGS  GIVES  THE  MINIMUM  WEIGHT  DES 
:1GN,  THIS  TABLE  MAY  ’/  5X.’SHOW  THAT  FROM  A  PRACTICAL  DESIGN  STAND 
:POINT  CONSIDERABLY  FEWER  RINGS  ’/  5X,’WILL  RESULT  IN  ONLY  A  SLIGHT 
;  WEIGHT  PENALTY.  ’) 

IF  (DRINGS.lt.  1.)  DRINGS  «  I. 

IF  (RINGSI.GT.O.)  RINGSI  «  RINGSI  4  DRINGS 

RINGS  «  RINGSI 

NRINGS  -  RINGS 

WRITE  (6,699)  ITERAT 

ITERAT  -  0 


C  . . .  NOTE!  . . . . . . . 

C  THE  NEXT  LINE  IS  USED  TO  FORCE  RETURN  TO  MAIN  PROGRAM. 

IF  (ITERAT.GE.0)  GO  TO  999 

C  OTHERWISE.  AN  INFINITE  LOOP  RESULTS  AND  FILE  -/THESIS  DAPSOUT  T1 
C  IS  FILLED  TO  C>\PACITY  WITHOUT  EVER  RETURNING  TO  MAIN  THE  PROGRAM. 
Q  j^oTE! . . . . 

IF  (RINGSI.LE.RINGSF)  GO  TO  3 
997  RINGSI  =  0. 

GO  TO  100 
999  CONTINUE 
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C  CLOSE  (UNrr»5,STATUS-’KEEP  ) 

C  CLOSE  (UNIT«6.STATUS»’KEEP) 

END 

SUBROUTINE  CUBlC(P,An.Bn^l2,B12,AI3,BI3,A22,B22,A23,B23, 

1  A33,B33) 

PI  ■  4.*ATAN(I.) 

El  •  -B11*B22*B33  +  B23**2*B11  +  B33*B12**2  -  2.*B12*BI3*B23 
I  +  BI3**2*B22 

FI  -  AI1*B22*B33  +  A22*B1I*B33  +  A33*BII*P22  -  2.*A23*B23*Bn 

1  -  AII*B23**2  -  2.*AI2*B12*B33  -  A33*B12**2  +  2.*(A12*B13*B23 

2  +  AI3*BI2*B23  +  A23*B12*BI3)  -  2.*AI3*B22*BI3  -  A22*B13**2 
G1  «  -A22*A33*B11  -  A11*A33*B22  -  AI1»A22*B33  +  A23**2*B11 

1  +  2.*A23*A11*B23  +  A12**2’B33  +  2.*AI2*A33'BI2  -  2.*(A12 

2  •A13*B23  +  AI2*A23*B13  +  A23*A13*B12)  +  A13**2*B22  +  2.*AI3 

3  •A22*B13 

HI  =  A11*A22*A33  -  A11*A23**2  -  A12**2*A33  +  2.*A12‘A13*A23 
1  -  A13**2’A22 
B1  «  Fl/El 
Cl  =  Gl/El 

R1  «  S0RT(ABS(C1  -  Bl‘*2/3.)**3/27.) 

THI  «  ACOS((Bl*Cl/3.-2.*Bl**3  y27..Hl/EI)/2./Rl) 

XI  =  2.*Rl**(l./3.)*COS(THl/3.)-Bl/3. 

X2  -  2.*Rl**(l./3.)*COS((THl  +  2.*PI)/3.)  -  Bl/3. 

X3  -  2.*Ri**(i./3.)*COS((THI  4.*PI)/3.)  -  BI/3. 

P  -  AMIN1(ABS(X1).ABS(X2)^S(X3)) 

RETURN 

END 

SUBROUTINE  CURVE(STRESS,STRAIN,NPNTS.Y,DYDX»E) 

C  THIS  SUBROUTINE  CALCULATES  THE  TANGENT  MODULUS  FOR  ANY 
C  VALUE  OF  STRESS  ON  A  CURVT  AS  REQUESTED  BY  SUBROUTINE  MODULI. 
DIMENSION  STRESS(NPNTS).STRAIN(NPNTS) 

COMMON  /COUNT/  ITERAT,INDIC,INDIC2,INDIC3, 

1  INDIC5.NDESIGN, ORTHO, NOKG,NCURVE,OK,ITERAM,OKM.ITERAG,OKG,ISEC, 

2  NPLOT,NOKM.NRINGS,NEGLAM,INDIC7 
ISAVE  «  NPNTS 

DO  100  I  «  1, NPNTS 
ITERAT  -  ITERAT  +  1 
IF  (Y.LE.STRESS(I))  ISAVE  =  I 
IF  (ISAVE.EO.I)  GO  TO  150 
100  CONTINUE 
150  CONTINUE 

IF  (ISAVE.OT.1)  GO  TO  200 
DYDX  -  E 
RETURN 
200  CONTINUE 

C  AT  THIS  POINT,  ISAVE  IS  AT  LEAST  2  AND  AT  MOST  10  . 

ERIGHT  -  (STRESS{ISAVE)  -  STRESS(ISAVE-1))  /  {STRAIN(ISAVE) 

1  -  STRAIN(ISAVE-l)) 

IF  (ISAVE.EO.NPNTS  AND.  Y.GE.STRESS(NPNTS))DYDX-AMAX1(ERIGHT,I0.) 
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IF  (ISAVE.EO.NPNTS  AND.  Y.GE.STRESS(NPNTS))  RETURN 
C  AT  THIS  POINT,  DYDX  HAS  BEEN  CALCULATED  AND  RETURNED  FOR  ANY 
C  STRESS  VALUE  LT.  THE  FIRST  POINT  OR  .GT.  THE  LAST  POINT  ON  THE 
C  INPUT  STRESS-STRAIN  CURVE. 

IF  (ISAVE.EQ.2)  ELEFT  =  E 

IF  (ISAVE.GT.2)  ELEFT  =  {BRIGHT  +  (STRESS(ISAVE-I)  -  STRESS(ISAVE. 

1  2))  /  (STRAIN(ISAVE-I)  •  STRAIN(ISAVE-2)))  /  2. 

IF  (ISAVE.LT.NPNTS)  ERIGHT  =  (BRIGHT  +  (STRESS(ISAVE+1)  - 
1  STRESS(ISAVE))  /  (STRAIN(ISAVE+I)  -  STRAIN(ISAVE)))  /  2. 

IF  (ERIGHT.GT.ELEFT)  GO  TO  300 

C  ELEFT  AND  ERIGHT  DEFINE  A  STRAIGHT  LINE  WITH  DYDX  ON  THIS  LINE. 
SLOPE  «  (ERIGHT  -  ELEFT)  /  (STRESS(ISAVE)  -  STRESS(ISAVE-l)) 

DYDX  =  ERIGHT  +  SLOPE*(Y  -  STRESS(ISAVE)) 

IF  (DYDX.LT.10.)  DYDX  =  10. 

RETURN 
300  CONTINUE 

WRITE  (6,60)  STRESS(1),STRAIN(1) 

60  FORMAT  ( ////  7X,’DAPS3  USER  FATAL  ERROR  MESSAGE  . .’  //  7X, 

1  ’THE  SLOPE  OF  THE  STRESS-STRAIN  CURVE  BEGINNING  *  //  7X, 

2  ’WITH  THE  POINT  (  ’,IPE10J,  ’ ,  *,E10.3,  *  )  IS  ’  //  7X. 

3  ’NOT  MONOTONICALLY  DECREASING  EVERYWHERE.  ’  ) 

STOP  1000 

END 

SUBROUTINE  GENERL(ETX,ETY,EF.VLT.VTL,PGEN,NSAVEP) 

C  THIS  SUBROUTINE  USES  THE  RENKEN  (RENZl-KENDRICK)  EQUATIONS  TO 
C  CALCULATE  THE  GENERAL  INSTABILITY  PRESSURE. 

C  FOR  THE  DERIVATION,  SEE  RENZI,  JOHN  R.,  "AXISYMMETRIC  STRESSES 
C  AND  DEFLECTIONS,  INTERBAY  BUCKLING,  AND  GENERAL  INSTABILITY  OF 
C  ORTHOTROPIC,  HYBRID,  RING-STIFFENED  CYLINDRICAL  SHELLS  UNDER 
C  EXTERNAL  HYDROSTATIC  PRESSURE,"  NAVAL  SURFACE  WEAPONS 
C  CENTER,NSWC 
C  TR  80-269,  19  APRIL  1981. 

C  EFE,  VXY,  AND  VYX  ARE  ALWAYS  ELASTIC  VALUES  HERE.  EF,  VLT,  AND 
C  VTI-  ARE  PLASTIC  VALUES  HERE  UNLESS  THIS  ROUTINE  IS  CALLED  DIRECTLY 
C  BY  RENKEN. 

REAL  L,IR,K1,K2,K3,KPR1ME,N,N2,L2,LF,LF2,LSA(5,5),B(5,5), 

1  BETA(5),WK(50).RALFA(10),RZ(50),RLAM(10),PLAM(11),VRZ(55),JAY 
COMMON  AR,IR,YR,LF,LS,ARM,L,RO,T,GLT,EX,EY,EFE,VXY,VYX,VPLAS,R 
COMMON  /COUNT/  ITERAT,INDIC,INDIC2,INDIC3, 

1  lNDlC5,NDESlGN,ORTHO,NOKG,NCURVE,OK,ITERAM,OKM,ITERAG,OKG,lSEC, 

2  NPLOT.NOKM,NRINGS,NEGLAM,INDIC7 

COMMON /GRAPH/ET(28),SIG(28),PB(28),PS(28),ES(28)AAA(6), 

1  BBB(6),CCC(8),DDD(6),EEE(5),FFF(7),GGG(4),HHH(7),000(8),ER(28), 

2  SIGF(28),SHAPEG(5),SHAPEI(3),WX(101) 

COMMON  /LAYER/  A1 1  A22,D1 1,D22 
COMPLEX  ALFA(5),Z(5,5),LAMBDA(5) 

EQUIVALENCE  (ALFA(1),RALFA(1)),(Z(I,I).RZ(1)),(LAMBDA(1),RIAM(1)) 
double  precision  da(.^,S),  db(5,5),  dbeta(S) 
complex*  16  dz(5,5),  dalfa(5) 

PI  =  4.*ATAN(1.) 
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PI2  -  PI‘*2 
PI3  -  PI**3 
R  «  RO  -  T/2. 

R2  -  R**2 
R3  -  R**3 
ALPHA  «  RO/R 
AL2  =  ALPHA**2 
T2  «  T*‘2 
T3  =  T**3 
ECC  =  YR  -  T/2. 

RN  •  NRINGS 

K1  s  4.*(RN  +  i.)**2/((2.*RN  +  3.)*(2.*RN  +  1.)) 

K2  •=  1.  +  T2/(12.*R2) 

K3  =  2. .  K2 
V7  «  VXY’VYX 
IF  (NDESIGN.EQ.-l)  GO  TO  200 
CYE  »  EY*T/(1.-V2) 

DXE  =  EX*T3/(12.»(1.-V2)) 

VT2  »  VLPVTL 
CX  «  ETX*T/(1.-VT2) 

CY  «  ETY*T/(1.-VT2) 

DX  -  ETX*T3/(12.*(1.-VT2)) 

DY  =  ETY*T3/(12.*(1.-VT2)) 

200  CONTINUE 

IF  (NDESIGN.NE.-1)  GO  TO  210 
CX  «  A22 
CYE=  All 
CY  »  All 
DXE  =  D22 
DX  «  D22 
DY  »  Dll 
210  CONTINUE 

C  PRINT  •,  •  GENERL  ".ETX.ETY.EF.VLT.VTL.GLT.T 
C  PRINT  •,  ’  ",EX,EY,EFE,VXY,VYX,RO,L 

C  PRINT  *, "  ".AR.IR,YR,LF,LS.RN 

L2  =  L*‘2 
LF2  «  LF**2 
BB  »  LF  -  LS 

C  FROM  HERE  THRU  GAMFE,  ELASTIC  PROPERTIES  MUST  ALWAYS  BE  USED. 
C  REMEMBER  THAT  V2  IS  ELASTIC. 

THETA  =  SORT(SORT(CYE»(l.-  V2)/(4.*DXE*R2)))*LS 

RR  »  RO  -  YR 

AEFF  =  AR*(R/RR)**2 

KPRIME  =  EFE/R2*(AEFF  +  BB*T) 

APRIME  *  KPRIME*R2/(CYE*(l..  V2))*(I.  •  VYX*AL2/2.)  - 
1  BB*(1.  -  VXY*ALPHA/2.) 

DPRIME  =  SINHfTHETA)  +  SIN(THETA) 

FPRIMl  «  2./(THETA*DPRIME)*(COSH{THETA)  -  COS(THETA)) 

GAMSE  »  CYE*(1.-  V2)»APRIME/(R*(KPR1ME/FPRIM1  + 

1  4.*THETA**4*DXE/LS*’3)) 
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BET  »  EFE*(AEFF  +  VYX*AL2/2.)/CYE/(l.-  V2)  - 

1  BB*(1.  -  VXY‘ALPHAy2.) 

GAMFE  «  -T/CYE/(I.-  V2)*((VYX*AL2/2.  -  1.)  + 

1  BET/(EFB*(AEFF  +  BB*T)/CYE/(1.- V2)  +  LS*FPRIM1)) 

EOR  =  ECC/R 
NODICE  •  0 

IF  (NSAVEP.EQ.-lOO  .OR.  (INDIC7.EQ.0  .AND.  NDESIGN.GT.O)) 

1  NODICE  «  1 
NSAVEP  »  2 
DO  81  K  =  2,12 
ITERAT  «  ITERAT  +  1 
N  =  K 
N2  =  N**2 

C  SETTING  JAY  »  0.  IS  THE  SAME  AS  ELIMINATING  THE  EIGHTH  INTEGRAL  IN 
C  EQUATION  (3-56)  ON  P.  3-23  OF  NSWC  TR  80-269  BY  RENZI.  KENDRICICS 
C  ISOTROPIC  SOLUTION  DID  NOT  HAVE  THIS  ANALOGOUS  INTEGRAL.  NOR  DID 
C  HE  HAVE  THE  EXACT  DEFINITIONS  FOR  GAMSE  AND  GAMFE  DEFINED  ABOVE 
C  AND  BY  EQUATIONS  (3-30)  AND  (3-44).  IN  ADDITION,  THIS  IS  THE  MOST 
C  EXACT  CLOSED  FORM  ORTHOTROPIC  SOLUTION  PRESENTLY  KNOWN. 

JAY  »  1. 

A(l,l)  »  CX*PI3*R/4./L  +  GL'PPI*'PL*N2*K2/4yR 

B(l,l)  «  PI‘L*N2/4.*(I.  -  GAMSE/R)  +  EF*PI*AR*G/vMFE‘(RN  - 

1  l.)*N2/4/r(l.  -»■  EOR) 

A(l,2)  =  -PI2*N/4.*(GLT*T*K3  +  VTL'CX) 

A(2,l)  «  A(l,2) 

B(l,2)  =  PI2*R»N/8. 

B(2,l)  »  B(l,2) 

A(l,3)  =  PI2/4.*(VTL*CX  -  GLT*T3*N2/67R2) 

A(3.1)  *  A(l,3) 

B(l,3)  «=  -Pl2*R/4. 

B(3,l)  «  B(I,3) 

A(2,2)  «  CY*PI*L*N2/4./R  +  GL'rPI3‘T*R*K2/4./L  -H 
I  EF*P!*AR*(RN  +  I.)*N2/4./R 

B(2,2)  *  PD‘R2*AL2/8./L  EF*PrAR*GAMFE*(RN  +  l.ymf2.mAY 
A(2,3)  =  -PI»N/4yR»(CY*L  -^  GLT*Pr2»T3/67L 
I  EF*AR*(RN  l.)*(l.  +  EOR*(I.  -  N2))) 

A(3,2)  =  A(2,3) 

B(2,3)  =  PrN/4./RyT*(L*TGAMSE  -» 

1  EFV-R*R*GAMFE*(RN  +  l.)*(JAY*(EOR*(N2-I.)-l.)-l,-EOR)) 

B(3,2)  =  B(2,3) 

A(33)  «  PI*L/i,/R*(CY  +  DY/R2*(N2  -  I.)**2  +  DX*PI2/L2*(PI2*R2/L2 

1  2.*VTL*(N2  -  1.)))  +  GLT‘PI3*T3*N2/12./R/L  +  EF*PI*AR*(RN  +  1.) 

2  /4./R*(l.+EOR»(EOR*(N2  -  l.)»*2  +  2.*(1.  -  N2)))  +  EF'PriR*(RN  + 

3  l.)*(N2  -  l.)**2/4./R3 

B(3,3)  «  PrL/4.»(EF*AR*GAMFE*N2T/L*(I.  EOR)*(RN  +  1.)  -  1. 

1+  PI2*R2*AJ^/2./U  N2/R*(R  -  GAMSE)) 

IF  (NDESIGN.EO.2.OR.NODICE.E0.1)  GO  TO  10 
C  ALWAYS  USES  THE  CUBIC  WHEN  NDESIGN  =  2. 

A(l,4)  =  -PI*N*K1*(GLT*T*K3  +  VTL’CX) 

A(4,l)  «  A(l,4) 
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B(l,4)  «  PI»R*N*Kl/2. 

B(4,l)  =  B(l,4) 

A(l^)  =  PrKl*(VYX*CX  -  GLT*  T3*N2/6yR2) 

A(5.1)  »  A(l^) 

B(l,5)  =  -PI’R’Kl 
B(5.1)  =  B(l.5) 

A(2.4)  =  CY*L*N2‘K1/R  +  GLrPI2*T*R*KI*K2/L 
A(4^)  =  A(2.4) 

B(2,4)  =  Kl*P12*R2*AL2/2./L 
B(4,2)  -  B(2.4) 

A(2^)  »  -CY*L*N*K1/R  - 
1  GLT*PI2*T3*N*KI/6yR/L 
A(5^)  »  A(2^) 

A(3,4)  =  A(2.5) 

A(43)  »  A(2^) 

B(2,5)  =  L*N*K1*GAMSE/R 
B(5,2)  =  B(2.5) 

B(3,4)  =  B(2,5) 

B(4,3)  =  B(2,5) 

A(3,5)  =  P12*Kl/LyR*(VTL*DX»(N2  -  I.)  +  DX»PI2*R2/L2  + 

1  GLT*T3*N2/3.)  +  K1*L/R*(CY  +  DY*(N2  -  l.)*((N2  -  l.)/R2 

2  VLT*PI2/L2)) 

A(5,3)  =  A(3^) 

B(3,5)  «  KrPI2*R2‘AL2/2./L  +  K1*L*(N2*(1.  -  GAMSETR)-!.) 

B(5,3)  »  B(3,5) 

A(4,4)  -  GLT*P13*T*R*L*K2/LF2  +  0.75*CY*PI*L*N2/R 

B(4.4)  »  PI3*L/2.*R2*AU/LF2 

A(4,5)  =  -PI*L*N/2./R*(1.5»CY  +  GL'rPI2*T3/3AF2) 

A(5,4)  =  A(4^5) 

B(4,5)  =  0.75*PI*L*N/rx*GAMSE 
B(S.4)  «  B(4.5) 

A(5,5)  =  PrL7./R*(1.5*CY  +  4.*DX*P12/LF2*(2.*PI2»R2/LF2  +  VTL*( 

1  N2  •  1.))  +  li*DY/R2»(N2  •  l.)**2)  +  GLT*PI3*T3*L*N2/3./R/LF2 
B(5.5)  =  PI3*L*R2*AL2/2./LF2  •  0.75‘PrL*(l.  +  N2*(GAMSE/R  -  1.)) 

C  DATA  LDA,LDB,IN,LDEVE05,5,5,5/  NEW  IMSL  NAMES 
DATA  lA.IB,IN.lZy5,5.5,5/ 

•  *  *  THE  IMSL  ROUTINE  GVCRG  (EIGZF)  IS  REQUIRED  HERE, 
do  8833  i  »  1,5 
DO  8833  J  =  1,  5 
DA(M)  =  A(I,J) 

DB(I,J)  «  B(U) 

8833  continue 

CALL  dGVCRG(IN,dA,IA,dB,IB,dALFA,dBETA,d4IZ) 

do  8834  i  <=  I,  5 
ALFA(I)  =  DALFA(I) 

BETA(I)  =  DBETA(I) 

DO  8834  J  =  1.  5 
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8834  z(i  j)  «  dz(ij) 


^  ••••••• 

C 

DO  48  I  -  1,IN 
ITERAT  =  ITERAT  +  1 
LAMBDA(I)  =  (0.,0.) 

48  IF  fBETA(I)  .NE.  0.)  LAMBDA(I)  =  ALFA(I)/BETA(I) 

C  PICK  OUT  THE  MINIMUM  EIGENVALUE  FOR  EACH  MODE  2  THRU  12  . 

ISAVc  e  1 

IF  (RLAM(1).LT.0.)  NEGLAM  «=  1 
IF  (RLAM(1).LT.0.)  RLAM(l)  »  500000. 

P  *=  RlAM(l) 

DO  80  I  *  1,9,2 
ITERAT  -  ITERAT  +  1 
IF  (RLAM(I).LT.O.)  NEGLAM  -  1 
IF  (RIAM(I).LT.O.)  RLAM(I)  »  500000. 

IF  (P.GT.RLAM(I))  ISAVE  =  I 
IF  (ISAVE.EQ.I)  P  *  RLAM(I) 

80  CONTINUE 
PLAM(K-l)  ^  P 

C  IN  A  BIG  VECTOR,  SAVE  ALL  THE  MODE  SHAPES  CORRESPONDING  TO  THE 
C  MINIMUM  EIGENVALUES  OF  EACH  MODE. 

DO  82  1  »  1,5 

ITERAT  ^  ITERAT  +  1 

VRZ(1+(K.2)*5)  a  RZ((ISAVE.1)‘5  +  1*2  -  1) 

82  CONTINUE 

IF  (NODICE.EQ.O)  GO  TO  49 

10  CALL  CUBIC  (P,A(1,1),P(1,1)A(1,2),B(1,2),A(13).B(13)A(2,2) 
1,B(2,2)A(2.3),B(2.3),A(33),B(33)) 

PLAM(K.l)  a  P 

49  CONTINUE 
IF(K.LT.3)  GO  TO  50 

IF  (P.LT.PLAM(K-2))  NSAVEP  =  K 
IF  (P.GT.PLAM(K  2))  GO  TO  83 

50  CONTINUE 

81  CONTINUE 

83  CONTINUE 

PGEN  =  PLAM(NSAVEP-.) 

IF  (NDESIGN.EQ.2.0R.NOD1CE.EO.I)  RETURN 
C  SAVE  THE  MODE  SHAPE  CORRESPONDING  TO  THE  ABSOLUTE  MINIMUM 
C  EIGENVALUE. 

DO  95  I  a  13 

ITFRAT  ■  ITERAT  +  1 

SHAPEG(I)  a  VRZ((NSAVEP  -  1)*5  -  5  +  I) 

95  CONTINUE 
RETURN 
END 

SUBROUTINE  IECAL(IE) 

C  THIS  SUBROUTINE  CALCULATES  THE  RING-SHELL  INERTIA  AS  CALLED  FROM 
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C  SUBROUTINES  RIBDW,  RIBAR,  ,\ND  RIBH.  (HERE,  AR  -  AREA.) 

COMMON  AR,IR.YR.LF,LS,DUM1(3),T.GXY.EX,EY,EF,VXY,VYX,VPLAS. 

1  DUM2(7),H,W,OD 

COMMON  /COUNT/  1TERAT,INDIC.INDIC2,INDIC3. 

1  INDIC5,NDESIGN, ORTHO, NOKG,NCURVE.OK,ITERAM,OKM,niERAG,OKG,ISEC. 

2  NPLOT,NOKM,NRlNGS,NEGLAM,INDIC7 
REAL  IR,1E,LF,LE,1S,LS 

AR  -  W  •  H 
IR  »  AR/12.*H**2 
VR  »  T  +  H/2. 

LS  =  LF  -  W 
RO  «  OD/2, 

R  -  RO  T/2. 

C  IF  (T.LE.0.)  PRINT  *,  *  lECAL  •,H.W,LF,LS,T,ITERAT,RO,R 
B  =  SQRT(SQRT(3.*(1..VXY*VYX)*EY/EX)/(R*T))*LF 
LE  -  (2.*(LF/B)*(COSH(B)-COS(B)))/(SINH(B)+SIN(B)) 

AS»LE*T 

IS-(AS/12.)*T‘*2 

YEAR  *  (T/2.*  AS+ YR*AR)/(AS+AR) 

RCG-RO.YBAR 

RR-RO-YR 

IE  «  EY/EF*(IS  +  AS*(R  •  RCG)**2)  +  IR  +  AR*(RR  •  RCG)**2 

RETURN 

END 

SUBROUTINE  ITERG(PIN,SECDEP.DW4>ASP) 

C  THIS  SUBROUTINE  PERFORMS  THE  DESIGN  ITERATIONS  ON  THE  RING 
C  HEIGHT,  WIDTH,  ASPECT  RATIO,  AND  SECIION  DEPTH  AS  REQUIRED. 
COMMON  AREA.IR,VR,LF,LS.DUM1(3),T,DUM2(8),P,DUM3(5),HEIGHT.WIDTH 
COMMON  /COUNT/  ITERAT.INDIC.INDIC2,INDIC3, 

1  INDICS,NDESIGN,ORTHO,NOKG,NCUBVE.OK,rTERAM,OKM,ITERAG,OKG,ISEC, 

2  NPLOT,NOKM,NRINGS,NEGLAM.INDIC7 
INTEGER  OKG 

REAL  IR.LF,LS 
ITERAT  -  ITERAT  +  I 

C  PRINT  *,"  ITERGl  ’,DUM I (3),TAREA,YR, HEIGHT, WIDTH, LS 
Q  ■  (P-PIN)/AMAXI(P.PIN) 

IF  (ABS(O).LT..001)  OKG  -  1 
C  PRINT  ITERGQ  ",P,PIN,Q,OKG 
IF  (OKG.EO.l)  RETURN 
ITERAG  ■  ITERAG  +  1 
HEIGHT  ■  HEIGHT/(1.  -  0/3.) 

C  IF  (T.LE.0.)  PRINT  *. '  ITERG  *,HEIGHT.WIDTH,LF,LS,T, ITERAT, Q, 

C  1  ITERAG 

IF  (ISEC.E0.1)  CALL  RIBH(SECDEP) 

C  PRINT  *:  RIBH  O  '.DUM I  (3),TAREA YR,HE1GHT.WIDTH,LS 

IF  (ISEC.EQ.0  .AND.  DW.GT.O.  .AND.  WIDTH.GT.DW)  CAIL  RIBDW(DW) 

AR  »  HEIGHT/WIDTH 
IF  (AR.LT.DASP)  CALL  RIBAR(DASP) 

C  PRINT  *,"  RIBARO  ",DUM1  (3), T,AREA.YR, HEIGHT, WIDTH, LS 
AR  -  HEIGHT/WIDTH 
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IF  (AR.GT.4.^D.WIDTH.LT..125)CALL  RIBDW(0.125) 

IF  (ITERAG.EQ.IOO)  OKG  «  1 
IF  (ITERAG.EQ.100)  NOKG  -  NOKG  +  1 
AREA  »  WIDTH  *  HEIGHT 
IR  -  AREA/I  2.*HEIGHT*  *2 
YR  -  T  +  HEIGHT/2. 

LS  »  LF  .  WIDTH 

C  PRINT  •;  ITERG2  ",DUM  1  (3),Tw^EA. YRJIEIGHT.WIDTHJ^ 

RETURN 

END 

SUBROUTINE  ITERM(PIN) 

C  THIS  SUBROUTINE  PERFORMS  THE  DESIGN  ITERATIONS  ON  THE  SKIN 
C  THICKNESS  AS  CALLED  FROM  THE  MAIN  PROGRAM. 

COMMON  AREA,IR,YR.DUM1  (4), RO,T,OUM2(7).R4>,DUMJ(5), HEIGHT, WIDTH 
COMMON  /COUNT/  ITERAT.INDIC,INDIC2,INDIC3. 

1  INDIC5,NDESIGN,ORTHO.NOKG,NCURVE.OK,ITERAM,OKM,ITERAG,OKG,ISEC, 

2  NPLOT.NOKM,NRINGS.NEGLAM.r.^IC7 
INTEGER  OKM 

REALIR 

C  PRINT  *,•  ITERM 1  ',RO,T.AREAY.R.HEIGHT,WIDTH,DUM  1(2) 

HERAT  •  ITERAT  +  1 
Q  -  (P-PIN)/AMAX1(P,P1N) 

IF  (ABS(0).LT..0008)  OKM  •  1 
C  PRINT  ITERMO  ".P.PIN.G.OrCM 
IF  (OKM.EO.I)  RETURN 
ITERAM  ■  ITERAM  +  1 
HI  -  T 

T  -  T/(l. .  Q/3.) 

R  -  RO  -  T/2. 

C  IF  (T.LE.0.)  PRINT  *,  *  ITERM  *, HEIGHT, WIDTH, DUMl(l),DUMi(2' 

C  1  T,ITERAT,Q.nERAM 

IF  (ITERAM.EQ.lOO)  OKM  «  I 

IF  (ITERAM.EQ.lOO)  NOKM  -  NOKM  +  1 

IF  (AREA.LE.0.)  RETURN 

HEIGHT  -  HEIGHT  +  HI  -  T 

AREA  ■  WIDTHMIEIGHT 

IR  -  AREA/I2.*HEIGHT**2 

YR  =  T  +  HEIGHT/2. 

C  PRINT  ITERM2  ",RO.T,AREA,YR,HEIGHT,WIDTH,DUMI(2) 

RETURN 

END 

SUBROUTINE  LOCAL(ETX,ETY,VXY,VYX,PREN,NSAVE) 

C  THIS  SUBROUTINE  USES  THE  RENVON  (RENZI-VON  MISES)  EQUATIONS  TO 
C  CALCULATE  THE  INTERBAY  OR  MONOCOQUE  SHELL  INSTABILITY  PRESSURE. 
C  THE  DERIVATION  OF  THE  EQUATIONS  CAN  BE  FOUND  IN  NSWC  TR  60-269. 

C  (SEE  THE  REFERENCE  IN  SUBROUTINE  GENERL.)  VXY  AND  VYX  ARE  PLASTIC 
C  VALUES  HERE  UNLESS  THIS  ROU  HNE  IS  CALLED  DIRECTLY  BY  RENVON. 
COMMON  DUM  l(4),LS,DUM2(2),RO.T,OXY,DUM3(6),R 
COMMON  /COUNT/  ITERAT, INDIC,INDIC2,INDIC3, 

1  INDlC5,NDESIGN,ORTHO,NOKO,NCURVE.OK,lTERAM,OKM,ITERAG,OKG,ISEC, 
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2  NPLOT,NOICM,NRINGS,NEGLAM,INDIC7 
CX)MMON  /GRAPH/ET(28),SIG(28),PB(28),PS(28),ES(28)A^A(6), 

1  BBB(6).0CC(8),DDD(6).EEE(5)4TF(7),GGG(4),HHH(7), 000(8), ER(28), 
2SIGF(28),SHAPEG(5),SHAPEI(3),WX(101) 

COMMON  /LAYER/  A11A22,D1 1,D22 
REAL  N,N2.N3.N44-S 
PI  ■  4.*ATAN(I.) 

PSAVE  -  1000000. 

R  =  RO  -  T/2, 

R2  =  R**2 

AA  -  PI'R/LS 

AA2  -  AA**2 

T3  «  'r*3 

V2  •  VXY*VYX 

IF  (NDESIGN.EQ.-1)  GO  TO  200 

CX  -  ETX*T/(I.-V2) 

CY  -  ETYT/(I.-V2) 

DX  -  ErX*T3/(12.*(l..V2)) 

DY  -  ETY*T3/(12.*(1.-V2)) 

200  CONTINUE 

IF  (NDESIGN.GT.-l)  GO  TO  210 
CX  -  A22 
CY  -  All 
DX  »  D22 
DY  -  Dll 
210  CONTINUE 
A3  •  VYX*CX*AA 
DO  80  NN  -  2,20 
ITERAT  -  ITERAT  +  1 
N  -  NN 
N2  -  N“2 
N3  -  N**3 
N4  -  N*»4 

A1  -  CX*AA2  ■¥  GXY*T*N2 

A2  -  -VYX'CX'AA'N  •  GXY’TAA’N 

A4  -  CX*AA*(AA2  -  VYX*N2) 

A5  -  CY»N3  •  VYX*CX*AA2*N  +  GXY/6yR2*AA2’T3*N  ■¥  DY/R2»N3 
A6  -  CY*(VXY*AA2  -  N2)  -  GXY/6JR2*AA2»T3*N2  -  DY/R2*N4 
A7  -  CX*AA*(VYX*(N2  -  I.)  -  AA2) 

A8  -  VYX*DX/R2*AA2*N  +  GXY/6yR2*AA2*T3»N  +  CY*N*(VXY'AA2  •  N2+1.) 
A9  -  CY*(N2-VXY*AAM,)-DX/R2*AA2*(AA2+VYX*N2)-CfXY/67R2*AA2*T3*N2 
AlO  -  R*(AA2/2.  N2  -  1.) 

PREN  -  ((A7»(A3*  A5.A2*  A6) + A8*(A1  *  A6.A3*  A4))/(A1  •  A5-A2*  A4)-A9)/A10 
IF  (PREN.GT.PSAVE)  GO  TO  81 
PSAVE  -  PREN 
NSAVE  -  NN 
DEN  -  APA3  -  A2*A4 
SHAPEI(l)  -  (A3*A5  -  A2*A6)/DEN 
SHAPEI(2)  ■  (ArA6  -  A3*A4)/DEN 
80  CONTINUE 
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81  CONTINUE 
PREN  -  PSAVE 
SHAPEI(3)  -  I. 

RETURN 

END 

SUBROUTINE  M  IDBAY(PSTR,SlGXM,SIGTM.SIGHOT) 

C  THIS  SUBROUTINE  USES  SOME  OF  THE  RENSALP  EQUATIONS  TO  CALCU- 
C  LATE  THE  MIDBAY  MIDPLANE  STRESSES  DUE  TO  A  PRESSURE  PSTR,  AS 
C  CALLED  FROM  REDUCE. 

REAL  IR.LF,LS,NETAI^ETA2.NI  J<2,L 

COMMON  AR.IR,YR,LF.LS,ARM,L,R0,T,GXY£X,EY^,VXY,VYX,VPLAS,R 
COMMON  /COUNT/  ITERAT,INDICINDIC2,INDIC3, 

1  INDIC5.NDESIGN,ORTHO.NOKG.NCURVE,OK.rrERAM,OKM,ITERAG,OKG,ISEC, 

2  NPLOT.NOKM.NRINGS.NEGLAM.INDIC7 
COMMON  /I  AYER/  AnA22,Dn,D22 

T3  -  T**3 
V2  -  VXY'VYX 
R2  -  R**2 

IF  (NDESIGN.EO.-l)  GO  TO  200 
CY  -  EY*T/(1..V2) 

DX  -  EX*T3/(12.*(1.-V2)) 

200  CONTINUE 

IF  (NDESIGN.GT.-l)  GO  TO  210 
CY  -  All 
DX  -  D22 
210  CONTINUE 
B-  LF-LS 
RR  -  RO  .  YR 
RATIO  ■  RO/R 

IF  (PSTR.LT.0.)  RATIO  •  (RO  T)/R 

THETA  -  LS*(CY*(1.-V2)/(4,*DX*R2))**.25 

GAMMA  -  PSTR*R2*RATIO**2/(4.*SQRT(DX*CY*(l.-V2))) 

IF  (OAMMA.GE.1.)  GAMMA  -  .999 
NETAl  -  .5*SORT(l,  -  G/VMMA) 

NETA2  -  .5*SORT(l.  +  G/VMMA) 

AEFF-AR*(R/RR)**2 
IF  (AR.LE.0.)  AEFF  ■  /VRM 

BETA  -  EF*(AEFF  +  B*T)»(l..VYX/2.*RATIO**2)/(CY*(l.-V2))  • 

1  B*(l.-VXY/2.»RATIO) 

SIGMAU  «  -PSTR*R/r 
N1  -  NETArTHETA 
N2  -  NETA2*THETA 

DEN  -  COSH(NI)*SINH(Nl)  /NETAl  -f  COS(N2)*SIN(N2)  /NETA2 

FI  -  (COSH(Nl)**2  -  COS(N2)**2)*4/rHETA/DEN 

F2  -  (COSH(Nl)*SIN(N2)/NETA2  SINH(Nl)*COS(N2)/NETAl)/DEN 

F4  ■  (■*-OOSH(Nl)*SIN(N2)/NETA2.SlNH(Nl)»COS(N2)/NETAl)/DEN 

DELTA  -  EF*(AEFF  +  B*T)  +  CY*(1.-V2)*LS*F1 

SIGXBM  -  6.*SIGMAU*BETA'F4/frDELTA)*SORT(DX*CY*(l..V2)) 

SIGTMF  -  SIOMAU*BETA/DELTA*CY*(l,-V2) 

SIGXO  -  .5*SlGMAU*RATIO*»2  +  SIGXBM 
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SIGXI  -  .3*SIGMAU*RATIO*»2  -  SIGXBM 
SIGTD  »  SIGMAU  -  SIGTMPF2  +  VYX'SIOXBM 
SIGH  »  SIGMAU  -  SIGTMF*F2  •  VYX'SIGXBM 
SIGXM  «  ABS((SIGXO  S10XI)/2.) 

SIGTM  ■  ABS((SIGTO  -f  S]GTI)/2.) 

QSTAR  -  .PSTR*B*(l..VXY/2.*RATIO)  +  SIGTMF*T*LS*F1/R 
SIGHOT  -  ABS(QSTAR*R/(AEFF  +  B*T)*(R/RR)) 

RETURN 

END 

SUBROUTINE  MODULI(X.YJ?.ETX.ETY,ETF.VLT.VTL,NCURVE) 

C  THIS  SUBROUTINE,  WITH  THE  HELP  OF  CURVE,  CALCULATES  THE  TANGENT 
C  MODULI  AS  REQUIRED  BY  SUBROUTINE  REDUCE.  THE  PLASTIC  POISSON 
C  RATIOS  ARE  ALSO  CALCULATED  HERE  FOR  USE  IN  THE  STABILITY  ROUTINES. 
COMMON  DUMI(I0),EX.EY,EF.VXY,VYXVPLAS,DUM2(4),SIGY 
COMMON  /EMOD/  STRESX(10),STRINX(I0),STRESY(I0).STRINY(I0), 
lSTRESF(10).SrRINF(IO),RATIOX,RATIOY,NXPNT,NYPNT,NFPNT 
ETX  -  EX 
ETY  -  EY 
ETF-  EF 
X  •  X/RATIOX 
YSAVE  ■  X 

IF  (X.LT.STRESX(1))  GO  TO  10 

CALL  CURVE(STRESX.yrRINXNXPNT,X,ETX,EX) 

10  X  -  X'RATIOX 
IF  (NCURVE.E0.1)  ETY  -  ETX 

IF  (NCURVE.EO.I  .AND.  F.GT.STRESX(1))  CALL  CURVE(STRESX 
1  STRINXNXPNT.F.ETF.EX) 

IF  (NCURVEEO.l)  GO  TO  30 
Y  •  Y/RATIOY 
YSAVE  •  Y 

IF  (Y.LT.STRESY(I))  GO  TO  20 

IF  (NCURVE.GT.1)  CALL  CURVE(STRESY,STRINY,NYPNT,Y,ETY,EY) 

20  Y  -  Y’RATIOY 

IF  (NCURVE.EQ.2  .AND.  F.OT.STRESY(l))  CALL  CURVE(STRESY, 

1  STRINY.NYPNT.F,ETF,EY) 

IF  (NCURVE.E0.2)  GO  TO  30 
IF  (F.LT.STRESF(1))  GO  TO  30 

IF  (NCURVE.E0.3)  CALL  CURVE(STRESF,STRINF,NFPNr,F,ETF,EF) 

30  CONTINUE 
IF  (ETX.OT.EX)  ETX  -  EX 
IF  (ETY.OT.EY)  ETY  -  EY 
IF  (ETF.OT.EF)  ETF  -  EF 

C  YI  IS  THE  PROPORTIONAL  LIMIT  STRESS  IN  THE  HOOP  DIRECTION. 

Y1  -  0. 

IF  (NCURVE.EO.I)  Yl  -  STRESX(I) 

IF  (NCURVE.GT.1)  Yl  -  STRESY(I) 

SLOPEV  -  0. 

IF  (Y1.LT.SIGY)  SLOPEV  -  (VPLAS  -  VYX)/(SIGY  -  Yl) 

B  -  VYX  •  SLOPEVYl 

C  LINEAR  VARIATION  OF  VTL  FROM  Yl  TO  SIGY. 
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VTL  -  SLOPEV*YSAVE  +  B 
IF  (YSAVE.lt. YI)  VTL  -  VYX 
C  BEYOND  YIELD,  VTL  IS  HELD  CONSTANT  AT  VPLAS. 

IF  (YSAVE.GT.SIGY)  VTL  »  VPLAS 
VLT  -  VTL*ETX/ETY 

C  WORDS  OF  WISDOM  ...  THE  PRODUCT  OF  A  STRESS  COMPONENT  AND  THE 
C  CORRESPONDING  STRAIN  COMPONENT  REPRESENTS  WORK  DONE  BY  THE 
STRESS. 

C  THE  SUM  OF  THE  WORK  DONE  BY  ALL  STRESS  COMPONENTS  MUSTBE  POSITIVE 
C  IN  ORDER  TO  AVOID  THE  CREATION  OF  ENERGY. ...  (P.  42,  JONES) 

IF  (VLT.GT.SQRT(ETX/ETY)/2.)  VLT  »  SORT(ETX/ETY)/l 
VTL  »  VLT*ETY/ETX 

IF  (VTL.GT.SQRT(ETY/ETX)/2.)  VTL  =  SORT(ETY/ETX)/l 
VLT  -  VTL*ETX/ETY 

C  PRINT  •,  •  MODULI  '.X.Y.F.ETX.ETY.ETF.VLT.VTL 
RETURN 
END 

SUBROUTINE  RCLOSE(RSAVE,WSAVE,RINGS, WEIGHT, RINGSI.DRINGS,RMAX. 

1  RCOUNT.OKR) 

C  THIS  SUBROUTINE  AIDS  IN  CLOSING  IN  ON  THE  OPTIMUM  NUMBER  OF 
C  RINGS  IN  THE  NDESiGN  -  2  OPTION. 

COMMON  AX)UNT/  nERAT,INDIC,INDIC2,INDIC3, 

1  INDIC5,NDESIGN,ORTHO,NOKG,NCURVE.OK,rTERAM,OKM,ITERAG,OKG,ISEC, 

2  NPLOT,NOKM,NRINGS,NEGLAM,INDIC7 
INTEGER  RCOUNT.OKR 

DIMENSION  RSAVE(RCOUNT).WSAVE(RCOUNT) 

IMIN  -  0 
RMIN  -  SOOOO. 

WMIN  «  SOOOO. 

DO  100  I  -  LRCOUNT 
ITERAT  •  ITERAT  +  1 
IF  (WSAVE(I)  .GE.  WMIN)  GO  TO  100 
IMIN  -  I 

RMIN  -  RSAVE(I) 

WMIN  >  WSAVE(I) 

100  CONTINUE 
RINGS  -  RMIN 
WEIGHT  -  WMIN 
IF  (DRINGS.LE.  1.)  GO  TO  300 
IF  (RMIN  .LE.  RINGSI)  GO  TO  200 
RINGS  -  RMIN  .  DRINGS 
RSAVE(l)  -  RINGS 
WEIGHT  -  WSAVE(IMIN-l) 

WSAVE(l)  -  WEIGHT 
200  RMAX  -  RINGS  +  2.*DRINGS  -  I. 

IF  (RMIN  .LE.  RINGSI)  RMAX  »  RINGS  +  DRINGS  -  1. 

DRINGS  ■  AINT(DRINGS/2.) 

IF  (DRINGS  .LT.  1.)  DRINGS  -  I. 

RCOUNT  -  1 
RETURN 
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300  OKR  -  2 
RETURN 
END 

SUBROUTINE  REDUCE(PINEL,NSAVE4) 

C  THIS  SUBROUTINE  CALCULATES  A  REDUCED  BUCKLING  PRESSURE  BASED  ON 
C  THE  STRESS  LEVEL  EXISTING  IN  THE  SHELL  AT  THE  SAME  PRESSURE. 
INTEGER  OK,OKM.OKG.ORTHO 
REAL  L,N,LF,IR.LS 

COMMON  ar.ir.yr,lf,ls.arm.l,ro,t,gxy.ex,ey,ef.vxy,vyx,vplas,r. 

1  PSTR.STR,SIGX.SIGY,SIGXl.SIGYl.DUM(2),OD 
COMMON  /COUNT/  ITERAT,INDIC.INDIC2,INDIC3, 

1  INDIC5,NDESIGN, ORTHO, NOKG,NCURVE,OK,ITERAM.OKM,ITERAG.OKG,ISEC, 

2  NPLOT,NOKM.NRINGS,NEGLAM,INDIC7 

COMMON  /GRAPH/ET(28).SIG(28).PB(28),PS(28),ES(28)AAA(6). 

1  BBB(6),CCC(8),DDD(6).EEE(5),FFF(7),GGG(4).HHH(7).000(8),ER(2S). 
2SIGF(28).SHAPEG(5),SHAPEI(3),WX(I0I) 

COMMON  /STRESS/  STRX,STRY,PIB,NSAVEI,PGEN,NSAVEG 
COMMON  /HENRY/  SRF.WM.WF 

PI  -  4.*ATAN(1.) 

INDIC6  =  0 

IF  (NCURVE  EQ.  4  .AND.  NPLOT  .EQ.  -2)  INDIC6  •  I 
C  PRINT  •,*  REDUCE  •.RO.T.AR.YR.DUM(I).DUM(2),LS 
C  IF  (INDIC6.EQ.1)  GO  TO  93 
PSI  •  0. 

SLOPES  «  0. 

SLOPEF  =  0. 

PBI  -  0. 

PB(1)  «=  0. 

IF  (J.EQ.0)  PEL  -  PINEL 
IF  (J.EQ.0)  J  =  I 

C  FOR  PROGRAM  OPTION  1  B.  RESORT  TO  THE  FOLLOWING  EFFECTIVE  MODULUS 
C  MUST  BE  USED.  (IN  OPTION  1 A  EX  -  E''  -  E.) 

E  -  SQRT(EX*EY) 

IF  (INDIC.EQ.1)  CALL  MIDBAY(PSTR,STRX,SrRY,SF) 

C  SF  HAS  AN  INFLUENCE  ON  THE  REDUCTION  OF  BOTH  ISOTROPIC  AND  ORTHO- 
C  TROPIC  RING  STIFFENED  SHELLS  SINCE  ETT  IS  USED  IN  THE  RING  TERMS. 

IF  (NDESIGN.GT.O  .AND.  OK.EQ.0)  INDIC5  -  1 
IF  (ORTHO.EO.O)  THEN 

STRX  «  SORT(STRX**2  +  STRV2  -  STRX'STRY) 

STRY  »  STRX 
END  IF 

IF  (ORTHO.EO.O)  THEN 
YIELDX  »  SIGY 
YIELDY  »  SIGY 
SIGX  »  SIGY 
END  IF 

IF  (ORTHO.EO.O)  GO  TO  47 

C  ALPHA  IS  A  MEASURE  OF  PERCENT  YIELD  CAUSED  BY  STRESSES  STRX  AND 
C  STRY.  YIELDX  AND  YIELDY  ARE  THE  STRESSES  WHICH  WILL  CAUSE  YIELD. 


107 


51  -  STRX/SIGX 

52  -  yrRY/SIGY 
ALPHA  -  AMAX1(SI.S2) 

YIELDX  -  STRX/ALPHA 
YIELDY  -  STRY/ALPHA 

47  CONTINUE 

C  (THE  NEXT  EQUATION  COULD  EQUALLY  BE  SLOPES  -  PSTR’YIELDX/STRX  .) 
IF(INDICEO.O.OR.  INDIC5.EQ.1 )  SLOPES  -  2.*PSTR/(STRX/YIELDX  + 

1  STRY/YIELDY) 

C  (INDIC  CANNOT  -  0  IF  INDICS  -  1  .) 

IF  (INDICS.EQ.l)  S1X)PEF  -  SF*YIELDX/STRX 
N  -  0. 

DS  -  0.04 
S  -  -0.04 
DO  90  I  - 

ITERAT  ■  ITERAT  +  1 
S  >  S  4-  DS 
SIG(I)  -  SMOO. 

C  THE  VALUES  AS  CALCULATED  IN  TT  '.P  NEXT  THREE  LINES  WILL  BE  RETAINED 
C  IF  (INDIC  .EQ.  0  .OR.  INDICS  .EQ.  I> . 

PS(I)  ■  SLOPES*S 
SX  -  S'YIELDX 
SY  -  S'YIELDY 
SF  ••  S*SLOPEF 

IF(INDICEQ.O  .OR.  INDICS.EQ.l)  GO  TO  9S7 
C  PI  IS  AN  ESTIMATE  OF  THE  PRESSURE  REQUIRED  TO  PRODUCE  S  PERCENT 
C  YIELD. 

PI  -  2.*T/R/SQRT(3.)*S»SQRT(YIELDX**2  +  YIELDY**2  -  YIELDX*YIELDY) 

IF  (S.EQ.0.)  GO  TO  956 

955  IF  (PI. LE  O.)  PI  ■  PB(I) 

ITERAT  ■  ITERAT  v  1 

CA  M!OU AY(P1,SX  OV.SF) 

XI  *  sy  vr>[SX'’*2  sv**2  -  sx*SY)/siax 

S)  •  .,.<  MGX 

S2»^'.s;jY 

IF(OnTHG.E0.1;Xl  -  AMAX1(S1.S2) 

ERO  »  (S-Xl)/AMAXI(S,X1) 

PI  -  Pl/(l.-ERO/2.) 

IF  (ABS(ERO).GT  .OOl)  GO  TO  955 

C  PI  IS  NOW  EXACTLY  (WITHIN  .1  PERCENT)  THE  PRESSURE  REQUIRED  TO 
C  PRODUCT  S  PERCENT  YIELD  AT  MIDBAY.  MID-PLANE  OF  THE  SHELL.  SF  IS 
C  THE  STRF.SS  IN  THE  FRAME  AT  PL 

956  CONTINUE 
PS(!)  -  PI 
SIGF(I)  -  SF 

957  CONTINUE 

IF  (ORTHO.EQ.O)  THEN 
SX  «  SORT(SX**2  SY*«2  -  SX'.SV^ 

SY  -  SX 
END  IF 
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CAU  MODULI(SX,SY,SF,ETX,ETY.ETF,VLT,VTUNCURVE) 

C  ARRAY  ET(I)  IS  AXIAL  MODULUS.  ARRAY  ES(I)  IS  HOOP  MODULUS. 

C  ER(I)  IS  RING  MODULUS. 

ET(I)  -  ETX 
ES(I)  -  ETY 
ER(I)  -  ETF 

C  THE  EFFECTIVE  MODULUS  IS  ONLY  USED  FOR  OPTION  I.B  . 

ETl  »  SORT(ETX*ETY) 

IF  (INDICEQ.0)  PB(I)  «  PEL’ETT/E 
INDICM  -  0 
INDIOG  -  0 
IF  (LEQ.l)  GO  TO  30 

IF  (ET(I).EO.ET(I-I)  .AND.  ES(I).EQ.ES(M))  INDICM  -  1 
IF  (INDICM.EQ.1  .AND.  ER(I).EO.ER(M))  INDICG  «  1 
30  CONTINUE 

IF  (INDICEQ.l  .AND.  J.EQ.l  .AND.  INDICM.EQ.O) 

I  CALL  LOCAL(ETX.ETY,VLT.VTL,PB(I),NSAVE) 

IF  (INDICEQ.l  .AND.  J.EQ.l  AND.  INDICM.EQ.1)  PB(I)  **  PB(I-l) 

IF  (J.EQ.2  .AND.  INDICG.EQ.O) 

1  CALL  GENERL(ETX,ETY,ETF.VLT.VTL,PB(I),NSAVE) 

IF  (J.EQ.2  .AND.  INDICG.EQ.I)  PB(I)  •  PB(I-I) 

IF  (I.EQ.28)  PB(I)  -  0. 

C  PRINT  *.  ■  REDUCEP  ‘.ETX,ETY.ETF,VLT.VTUPS(I),PB(I),NSAVE4 
IF  (N.EQ.1..0R.PB(I).GT.PS(I))  GO  TO  90 
N  -  1. 

SSAVE  -  SIG(M)/I00. 

IF  ('ND1C5.E0.1)  GO  TO  95 
90  CONTINUE 
95  CONTINUE 
IF  (INDIC.EQ.0)  THEN 
ETSAVE  «  E 
ESSAVE  -  E 
END  IF 
DS  «  0.004 
S  -  SSAVE  -  2.*DS 
N  -  0. 

DO  89  1-1,13 
ITERAT  -  ITERAT  +  1 
IF  (N.EO.L)  GO  TO  86 
S  -  S  4^  DS 
PSI  ■  SLOPES'S 
SX  -  S'YIELDX 
SY  -  S' YIELD  Y 
SF  -  S'SLOPEF 

IF  (INDIC.EO.O  .OR.  INDIC5.E0.1)  GO  TO  857 

PI  ■  2,*T/RySORT(3.)'S'SORT(YIELDX"2  +  YIELDY"2  -  YIELDX'YIELDY) 
855  IF  (PI  LE  D.)  PI  -  PBl 
ITERAT  -  ITERAT  +  1 
CALL  MIDBAY(PI.SX,SY,SF) 

XI  -  SORT(SX"2  +  SY"2  -  SX'SY)/SIGX 
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51  -  SX/SIGX 

52  »  SY/SIGY 

IF  (ORTHO.EQ.l)  XI  -  AMAX1(S1.S2) 

ERO  =  (S-X1)/AMAX1(S,X1) 

PI  »  Pl/(1.  -  ERO/2.) 

IF  (ABS(ERO).GT..001)  GO  TO  855 
PSI  B  PI 
857  CONTINUE 

IF  (ORTHO.EQ.O)  THEN 
SX  -  SQRT(SX**2  +  SY**2  -  SX'SY) 

SY  -  SX 
END  IF 

C  PRINT  •.  •  CALL  MODULI  FROM  REDUCE  NEAR  LABEL  85r 
CALL  MODULI(SX.SY,SF.ETX,ETY,ETF,VLT,VTL,NCURVE) 

C  THE  EFFECnVE  MODULUS  IS  ONLY  USED  FOR  OPTION  LB  . 

ETI  «  SQRT(ETX*ETY) 

ETSAVE  -  ETX 
ESSAVE  «  ETY 
EFSAVE  =  ETF 

IF  (INDICEQ.0)  PBI  «  PEL*ETI/E 
INDICM  «  0 
INDICG  «  0 
IF  (l.EQ.l)  GO  TO  40 

IF  (ETX.EQ.ETSAVE  .AND,  ETY.EQ.ESSAVE)  INDICM  -  I 
IF  (INDICM.EQ.1  AND.  ETF.EO  EFSAVE)  INDICG  «  1 
40  CONTINUE 

IF  (INDICEQ.l  .AND.  J.EQ.l  .AND.  INDICM.EQ.O) 

1  CALL  LOCAL(ETX,ETY,VLT,VTL.PBI,NSAVE) 

IF  (J.EQ.2  .AND.  INDICG.EQ.O) 

1  CALL  GENERL(ETX,ETY,ETF.VLT,VTL,PBI,NSAVE) 

IF  (I.EQ.13)  PBI  »  0. 

IF  (N.EQ.1..0R.PBI.GT.PSI)  GO  TO  87 
N  »  1. 

SLOPEB  »  n»BI  -  PSAVE)/DS 

SSLOPE  -  (PSI  -  SSAVE)/DS 

BINTER  -  PBI  -  SLOPEB'S 

SINTER  »  PSI  -  SSLOPE'S 

SIGMA  -  (SINTER  -  BINTER)/(SLOPEB  -  SSLOPE) 

PINEL  -  SSLOPE-SIGMA  +  SINTER 
87  PSAVE  «  PBI 
SSAVE  «  PSI 
89  CONTINUE 
86  CONTINUE 

C  IF  (T.LE.0.)  PRINT  *, "  REDUCE  ",T,ITERAT,PINEL,ITERAM,ITERAG 
IF  (INDIC5.EQ.1)  GO  TO  99 

IF  (INDIC.EO.l.AND.AR.GT.0,.ANDJ.EO.I)WRITE  (6.692)  PIB.NSAVEI 
IF  (INDIC,EO.l.AND.AR.LE.0„ANDJ.EQ.I)WRITE  (6,693)  PIB.NSAVEI 
IF  (J.E0.2)  WRITE  (6,655)  PGEN.NSAVEG 
SIGMA  -  SIGMA*  100. 

WRITE  (6,698)  PINEUNSAVE, SIGMA 
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IF  (INDICEQ.l  .AND.  J.EQ.l)  WRITE  (6,694)  SHAPEI 

IF  (J.E0.2)  WRITE  (6,662)  SHAPEG 

IF  (J.EQ.2  .AND.  NEGLAM.EQ.l)  WRITE  (6,663) 

C  •••••  PLOT  OUTPUT  REMOVED*  •••  * 

IF  (NPLOT.GT.O)  GO  TO  99 

655  FORMAT  (  /  5X,’ELACT1C  GENERAL  INCTAPILITY  PRESSURE  «  ’,IPE10J, 
15X,’MODE  =  ’,12) 

662  FORMAT  (  5X,’A1**,1PE10J,’  B1»*,EI0J,’  C1«’,E103,’  B2=’, 

1  E103,’  C2-’,E103) 

663  FORMAT  (  5X,’(A  NEGATIVE  EIGENVALUE  WAS  IGNORED.)’ ) 

692  FORMAT(  /5X,’ELAST1C  INTER-BAY  INSTABILITY  PRESSURE  -  ’,1PE103, 
13X,’MODE  «  ’,12) 

693  FORMAT  ( /SX,’ELASTIC  MONOCOQUE  SHELL  INSTABILITY  PRESSURE  »  ’, 

1  1PE1033X,’MODE  »  ’,12) 

694  FORMAT  (5X,‘BUCKLED  SHAPE  ...  •,2X,’A  »  ’,1PE103,2X,’  B  «  ’, 

1  E103,2X,’  C  »  ’,0PF5  J) 

698  FORMAT  (  5X,’***  INELASTIC  BUCKLING  PRESSURE  »  ’,1PE103,1X,’PSI 
l.’,4X,’MODE  =  ’,12/  5X,’THE  STRESS  LEVEL  IN  THE  SHELL  »  ’, 

2  0PF10.4,’  PERCENT  YIELD.’) 

801  FORMAT(’  INELASTIC  BUCKLING  PRESSURE  -  ’.IPEIOJ,’  PSI,  MODE  -  ’, 

1  12) 

802  FORMATC  THE  REPRESENTATIVE  STRESS  LFVEL  -  ’,ri0.4,’  PERCENT  YIELD 
1’) 

803  FORMATC  YIELD  STRESSES  LC  «  ’,F7.0,’,  LT  -  ’^.0,’,  HC  -  ’, 

1  F7.0,’,  HT  «  ’,F7.0,’  PSI  ’) 

804  FORMATC  SHELL  LENGTH  -  ’,F83.’  IN,  OD  •  T  -  ’4?53) 

805  FORMATC  ELASTIC  BUCKLING  PRESSURE  -  ’,1PE103,’  PSP) 

806  FORMATC  MIDBAY,  MIDPLANE  STRESSES  «  ’,1PE103,’  (HOOP)  AND  ’, 

1  ElOJ,’  (AXIAL)  ’) 

807  FORMATC  UNDER  A  PRESSURE  OF  ’,IPE103,’  PSI.  ) 

808  FORMAT(  ’YIELD  STRESSES  LC  -  ’ J7A’,  LT  =  ’,F7.0,’,  HC  *  ’, 

1  F7.0,’,  HT  »  ’,F7.0,’  PSI  ’) 

99  CONTINUE 
RETURN 
END 

SUBROUTINE  RENKEN(PINELG) 

C  THIS  SUBROUTINE  DIRECTLY  CALLS  GENERL  TO  CALCULATE  THE  ELASTIC 
C  GENERAL  INSTABILITY  PRESSURE  WHEN  INDIC5  =  0  .  WHFN  NCURVE  .NE.  4 
C  REDUCE  IS  CALLED  WHICH  CALLS  GENERL  TO  CALCULATE  THE  INELASTIC 
C  VALUE  BY  USING  TANGENT  MODULI. 

COMMON  DUM(10),EX,EY,EF,VXY,VYX 
COMMON  /COUNT/  nrERAT,INDIC,INDIC2,INDIC3, 

1  INDIC5,NDESIGN,ORTHO,NOKG,NCURVE,OK,rrERAM,OKM,ITERAG,OKG,ISEC, 

2  NPLOT,NOKM,NRINGS,NEGLAM,INDIC7 

COMMON  /STRESS/  STRX.STRY,P?B,NSAVEI,PGEN,NSAVEG 
J  «  2 

IF  (INDIC5.EO.0)  CALL  GENERL(EX,EY,EF,VXY,VYX,PGEN,NSAVEG; 

IF  (NCURVE.EQ.4  .AND.  NPLOT.E0..2)  CALL  REDUCE(PINEL,NSAVE3) 

C  PINEL  AND  NSAVE  IN  THE  ABOVE  CALL  IS  DUMMY  SINCE  PGEN  AND  NSAVEG 
C  IS  IN  COMMON  WITH  REDUCE, 
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IF  (NCURVE.EQ.4)  RETURN 
C  PRINT  *,  "  CALL  REDUCE  FROM  RENKEN* 

CALL  REDUCE(PINEL,NSAVE4) 

PINELG  «  PINEL 

RETURN 

END 

SUBROUTINE  RENSALP(PCOKT,SIGVOF,SIGVIF) 

C  THIS  SUBROUTINE  USES  NONLINEAR  STRESS  EOUA^ONS  (WHICH  INCLUDE 
C  THE  BEAM-COLUMN  NONLINEARITY  PRODUCED  BY  AXIAL  COMPRESSION  WITH 
C  AJCIAL  BENDING)  TO  CALCULATE  THE  AXIS YMMETRIC  ELASTIC  DEFORMATIONS 
C  AND  STRESSES  IN  A  RINO-SHFFENED  CIRCULAR  CYUKDRICAL  SHELL  UNDER 
C  HYDROSTATIC  PRESSURE.  THESE  SAME  EQUATIONS  ARE  USED  TO  CALCULATE 
C  THE  AXISYMMETRIC  COLLAPSE  PRESSURE.  THE  RENSALP  (RENZI  SALERNO- 
C  PULOS)  EQUATIONS  ARE  GOOD  FOR  ORTHOTROPIC  HYBRID  SHELLS  AS  WELL 
C  AS  ISOTROPIC  SHELLS. 

C  FOR  THE  DERIVATION,  SEE  RENZI,  JOHN  R..  ’OPTIMIZATION  OF 
C  ORTHOTROPIC,  NON-LINEAR,  RING-SnFFENED  CYLINDRICAL  SHELLS  UNDER 
C  EXTERNAL  HYDROSTATIC  PRESSURE  AS  APPLIED  TO  MMC  MATERIALS",  NSWC 
C  TR  79-305,  30  SEPTEMBER  1979.  SEE  ALSO  RENZI,  NSWC  TR  80-269. 

INTEGER  OKT,ORTHO 

REAL  IR.LF,LS,NETA1  ,NETA2.N  1,N2.L,KSI(I0) 

COMMON  AR.IR,YR,LF,LS.ARM,L,RO,T,GXY,EX,EY,EF,VXY,VYX,VPLAS,R, 

:  PSAV,SrR,SIGX,SIGY,SIGXl.SIGYl 
COMMON  /COUNT/  rrERAT,INDIC,INDIC2,lNDIC3, 

:  INDIC5,NDESIGN,ORTHO,NOKG,NCURVE,OK,rrERAM,OKM,lTERAG,OKG,lSEC, 

;  NPLOT,NOKM,NRINGS,NEGLAM,INDIC7 

COMMON  /EMOD/  STRESX(10),STRINX(10),STRESY(10),STRINY(10), 

:  STRESF(10),STRINF(IO),RATIOX,RATIOY,NXPNT,NYPNT,NFPNT 
COMMON  /LAYER/  A1 1,A22,D1 1,D22 
COMMON  /HENRY/  SRF,WM,WF 

C  WHEN  OKT  «  1,  RENSALP  IS  EXPRESSLY  USED  TO  QUICKLY  CALCULATE  THE 
C  AXISYMMETRIC  COLLAPSE  PRESSURE. 

C  PRINT  RENSALP  •,RO,T,AR,YR,LS 
T3  =  T”3 
V2  =  VXY’VYX 
IF  (NDESIGN.EQ.-I)  GO  TO  470 
CY  =  EY*T/(1..V2) 

DX  «  EX’T3/(12.’(1.-V2)) 

IF  (NDESIGN.GT.-l)  GO  TO  480 
470  CONTINUE 
CY  »  All 
DX  =  D22 
480  CONTINUE 
PSTR  =  PSAV 
B  e  LF  -  LS 
RR  -  RO  -  YR 
R  »  RO  -  T/2. 

R2  «  R”2 
RATIO  =  RO/R 
IND1C2  -  0 
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IF  (OKT.EQ.l)  GO  TO  490 
IF  (PSAV.LT.O.)  RATIO  -  (RO-T)/R 
IF  (NDESIGN.lt.  1)  WRITE  (6.60)  PSAV 
IF  (NDESIGN.EQ.l)  WRITE  (6,61)  PSAV 
490  CONTINUE 

THETA  -  LS*(CY*(1.-V2)/(4.*DX*R2))‘*.25 

AEFF»AR*(R/RR)**2 

IF  (AR.LE.0.)  AEFF  =  ARM 

BETA  «  EF*(AEFF  +  B*T)*(l.-VYX/2.*RATIO»*?.)/(CY*(l.-V2))  - 
1  B*(l.-VXY/2.*RATIO) 

600  GAMMA  »  PSTR*R2*RATIO**2/(4.-SQRT(DX*CY*(l..V2))) 

ITERAT  »  ITERAT  +  1 
IF  (GAMMA.GE.1.)  GAMMA  -  .999 
NETAl  =  .5’SQRT(1.  -  GAMMA) 

NETA2  -  .5*SQRT(1.  +  GAMMA) 

SIGMAU  «  -PSTR*R/r 
N1  «  NETArTHEFA 
N2  «=  NETA2*THETA 

DEN  »  COSH(Nl)*SINH(Nl)  /T-IETAl  +  COS(N2)*SIN(N2)  /NETA2 
FI  =  (COSH(Nl)**2-COS(N2)**2)*4yiTIETA/DEN 
F2  «  (COSH(Nl)*SIN(N2)/NETA2  +  SINH(Nl)*COS(N2)/NETAl)yDEN 
F4  «=  (+COSH(Nl)*SIN(N2)/NETA2  -  SINH(Nl)*COS(N2)/NETAl)/DEN 
DELTA  «  EF*(AEFF  +  B*T)  +  CY*(1.-V2)*LS*F! 

IF  (OKT.EQ.l  .OR.  INDia.EQ.l)  GO  TO  700 
F3  »  (.COSH(Nl)*SINH(Nl)/NETAl  +  COS(N2)*SIN(N2)/NETA2)/DEN 
F5  -  (COSH(Nl)»SIN(N2)/NETA!  •  SINH(Nl)*OOS(N2)/NETA2)/DEN 
FAC4  =  F2*COSH(Nl)*COS(N2)  +  F5*S1NH(N1)*SIN(N2) 

SIGXBM  «=  6.*S1GMAU*BETA‘F4/(T*DELTA)*S0RT(DX'CY»(1..V2)) 
SIGTMF  «  SIGMAU*BETA/DELTA*CY*(L-V2) 

FAC2  «  6.*SIG7MF*F3*SQRT(DX/(CY*(1,-V2)))A' 

SIGXO  =  .5»SIOMAU*RATIO**2  +  SIGXBM 
SIGXI  »  .3'SIGMAU*RATIO«*2  -  SIGXBM 
SIGTO  »  SIGMAU  •  SIGTMF*  F2  +  VYX'SIGXBM 
SIGTT  «  SIGMAU  -  SIGTMF*F2  -  VYX'SIGXBM 
SIGXM  »  (SIGXO+SIGXI)/2. 

SIGTM  »  (SIGTO +SIGTI)/2. 

SXI  «  ABS(SIGXI/S1GX) 

SXO  -  ABS(SIGXO/SIGX) 

Sn  «  ABS(SIGTI/SIGY) 

STO  »  A3S(SIGT0/SIG'0 
IF  (SIGXI.GT.0.)  SXI  ^  SIGXI/SIGXl 
IF  (SIGXO.GT.O.)  SXO  «  SIGXO/SIGXl 
IF  (SIGTI.GT.0.)  sn  »  SICn/SIGYl 
IF  (SIGTO.GT.0.)  STO  »  SlGTO/SIGYl 

C  KEEP  IN  MIND  THAT  FOR  ORTHO  *  0,  SIGX  -  SIGY  «  SlGXl  »  SlGYl. 
IF  (ORniO.EQ.O)  SIGVO«SQRT(SIGTO**2+SlGXO**2-SiGTO*SIGXO)/ 
1S1GX*100. 

IF  (ORTHO.EQ.O)  SIGVI  «  SORT(SIGn**2+SIGXl»-2-SlGn*SIGXI)  / 
ISIGX’IOO. 

IF(ORTHO.EQ.l)  SIGVO  »  AMAX1(SXO,STO)*100. 


IF(ORTHO.EQ.l)  SIGVI  »  AMAXl(SXI.Sn)*100. 

SIOXOF  =  .5*SlGMAU*RATIO**2  +  FAC2 
SIGXIF  -  .5*S:g:<AU*RATIO*‘2  -  FAC2 
SIGTOF  «  SIGMAU  -  SlGTMF  +  VYX*FAC2 
SIGTIF  SIGMAU  -  SlGTMF  -  VYX*FAC2 
SX:  =  ABS(SIGXIF/SIGX) 

SXO  =  ABS(SIGXOF/SlGX) 
sn  «=  ABS(SIGTIF/SIGY) 

STO  «  ABS(SIGTOF/SIGY) 

IF  (SIGXIF.GT.O.)  SXI  »  SIGXIF/SIGXl 
IF  (SIGXOF.GT.O.)  SXO  »  SlGXOF/SIGXl 
IF  (SIGTIF.GT.O.)  STI  »  SIOTIF/SIGYI 
IF  (SIGTOF.GT.O.)  STO  «  SIGTOF/SIGYl 

IF  (ORTHO.EQ.O)  SIGVOF  =  SQRT(SIGTOF**2+S1GXOF**2  SIGTOF*SIGXOF)  / 
ISIGXMOO. 

IF  (ORTHO.EO.O)  SIGVIF  «  SQRT(S1GTIF••2-^SIGXIF••^SIG^^F•SIGXIF)  / 
ISIGXMOO. 

IF(ORTHO.EQ.l)  SIGVOF  «=  AMAX1(SXO,STO)*100. 

IF(ORTHO.EQ.l)  SIGVIF  -  AMAXl(SXl.Sn)' 100. 

OSTAR  -  -PSTR*B*(l.-VXY/2.‘RATIO)  +  SIGTMF*T*LS‘F1/R 
SIGHOT  -  QSTAR‘R/(AEFF  +  B*T)*(R/RR) 

SRF  -  SIGHOT 

WINF  «  *PSTR*R2‘(!.  -  VYX/2.»RAT10**2)/(CY*(1..V2)) 

WM  »  WINF  +  PSTR*R2*BETA*F2/DELTA 
WF  -  WINF  +  PSTR*R2‘BETA*FAC4/DELTA 
SXI  »  ABS(SIGXM/SIGX) 

Sn  -  ABS(SIGTM/SIGY) 

IF  (SIGXM.GT.0.)  SXI  SIGXH^IGXl 
IF  (SIGTM.GT.0.)  STI  =  SIGTM/SIGYI 
IF  (ORTHO.EQ.O)  STR  -  SQRT(SIGXM-2+SIGTM*  '2 
1  •  SIGXM*SIGTM)/SIGX*100. 

!F(ORTHO.EQ.l)  STR  *  AMAX1(SXI,STI)*100. 

WRnE(6,100) 

WRITF(6,200)  SIGXO,SIOXI.SlGTO,SIGTI,SIGVO,SlGVI,WM 

WRITE(6,50)STR 

WRrrE(6,250) 

WRITE(6,200)  oIGXOF, SIGXIF, SIGTOF.SIGTIF, SIGVOF, S'GVIF.WF 

WRrrE(6,300)SIGH01 

IF  (ORTHO.EQ.O)  WRITE  ^6,320) 

IF  (ORTH  ).EQ.l)  WRITE  (6,322) 

WRITE  (6325)  S1GY,SIGX,SIGY1,SIGXI 
IF  (NCURVE.E0.4)  GO  TO  80 
WRITE  (6310) 

DO70I-l,NXPNT 

70  KSI(I)  «  STF.ESX(I)/1000. 

WRITE  (6311)  (KSI(I),1»1,NXPNT) 

WRITE  (6314)  (STRINX(I),1«1,NXPNT) 

IF  (NCURVE  EQ.l)  GO  TO  80 

DO  71  I  «  1,NYPNT 

71  KSI(I)  «  STRESY(1)/1000. 


A 
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WRITE  (6312)  (KSI(I).I-1,NYPNT) 

WRITE  (6314)  (STRINY(I),I»1,NYPNT) 

IF  (NCURVE.E0.2)  GO  TO  80 
DO  72  1  -  1,NFPNT 
72  KSI(I)  »  STRESF(I)/1000. 

WRITE  (6313)  (KSI(I),I-1,NFPNT) 

WRITE  (6314)  (STRINF(1).I-1,NFPNT) 

80  CONTINUE 
INDIC2  »  1 

IF  (PSAV.LE.0.)  RETURN 

700  A  -  6.*BETA*F4/(T*DELTA)*SQRT(DX*CY*(I.-V2)) 

R5  »  -1.  +  BETA/DELTA*CV(1..V2)*F2  •  VYX»A 
R6  »  -03*RATIO‘*2  -  A 

IF  (ORTHO.EQ.O)  PC  -  SIGX*T/R/SQRT(R5**2+R6**2-R5*R6) 

IF  (ORTHO.EQ.O)  GO  TO  701 
PCX  -  ABS(SIGX*T/R/R6) 

PCY  =  ABS(SIGY*T/R/R5) 

PC  -  AM1N1(PCX.PCY) 

701  CONTINUE 

ER  -  (PSTR-PC)/AMAX1(PSTR,PQ 
PSTR-PO(l.-ER/2.) 

IF(ABS(ER).GT..001)  GO  TO  600 
IF  (OKT.EQ.l)  RETURN 
WRITE  (6350)  PC 
RETURN 

50  FORMAT(  /3X,'CALCULATED  MIDBAY  MEMBRANE  STRESS  LEVEL  IS  ’,F10.4,* 

1  PERCENT  YIELD/) 

60  FORMAT  ( /SX’ANALYSIS  PRESSURE  IS  '3^.0/  PSI/) 

61  FORMAT  (  /5X,'DESIGN  PRESSURE  IS  *3^6.0,’  PSI/) 

100  FORMAT(  /5X/SHELL  STRESSES  AT  MIDBAV) 

200  FORMAT(  35X/OUTSIDE',12X,’1NSIDE7,10X,’LONGITUDINAL(PSI)’,4X, 
lF10.0,8X,F10.0/,10X,’aRCUMFERENTIAL  (PSI)MX,  F10.0,8X,FI0.0/, 

210X/PERCENT  YIELD  ••,6X3^10.4,8X3^10.4,  /,10X, 

3’RADIAL  DEFLECTION  (IN.)',6X,F7.4) 

250  ;  ORMAT(  /5X.’SHELL  STRESSES  AT  THE  FRAME  ) 

300  FORMAT(10X,’HOOP  STRESS  IN  FRAME  (PS1)’,F10.0) 

310  FORMAT  (  5X, ’STRESS-STRAIN  DATA  ECHO  ...  (STRESS  IS  IN  KSI  FOR  EC 
IHO  ONLY.)’) 

311  FORMAT  ( 5X.’X',10(1X,F6.2)) 

312  FORMAT  ( 5X,’r.l0(lXF6.2)) 

313  FORMAT  ( 5X/r,10(IXF63)) 

314  FORMAT  ( 6X,10(1X,F63)) 

320  FORMAT  (/5X,’*  YIELD  STRESSES  (USED  WITH  VON  MISES  YIELD  CRITERl 
ION)  ’) 

322  FORMAT  (/5X,’*  YIELD  STRESSES  (USED  WITH  MAXIMUM  S1RESS  YIELD  CR 
IITERION)  ’) 

325  FORMAT  (  24X.’HOOP  (Q  -  ’,F7.0,’  PSr,8X,’ AXIAL  (C)  -  ’,F7.0, 

1  ’  PSI  ’/24X,’HOOP  (T)  -  ’,F7.0,’  PS1’,8X,’AXIAL  (T)  -  ’,F7.0, 

2  ’  PSI  ’) 

350  FORMAT(  /5X,’AXISYMMETRIC  COLLAPSE  PRESSURE  «  ’.FO.O  PSI,  BASED 
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1  ON  OUTER  FIBER  V6X,’STRESSES  AT  MID-BAY.’) 

END 

SUBROUTINE  RENVON(PINELB) 

C  THIS  SUBROUTINE  DIRECTLY  CALLS  LOCAL  TO  CALCULATE  THE  ELASTIC 
C  INTERBAY  OR  MONOCOQUE  SHELL  INSTABILITY  PRESSURE  WHEN  INDICS  »  0 

C  WHEN  NCURVE  NE.  4,  REDUCE  IS  CALLED  WHICH  CALLS  LOCAL  TO  CALCU- 
C  LATE  THE  INELASTIC  VALUE  BY  USING  TANGENT  MODULI. 

COMMON  DUM(10),EX,EY,EF,VXY,VYX 
COMMON  /COUNT/  ITERAT,INDIC,INDIC2,INDIC3. 

1  INDICS, ndesign,ortho.nokg,ncurve,ok,iteram,okm.iterag,okg.isec, 

2  NPL0T,N0KM.NRINGS,NEGLAM,INDIC7 

common  /stress/  strx.stry,pib.nsavei.pgen.nsaveg 

IF  (INDIC5.EQ.0)  CALL  LOCAL(EX,EY,VXY,VYX,PIB,NSAVEI) 

IF  (NCURVE.E0.4)  RETURN 
J  »  I 

C  PRINT  *.  "  C/UJ-  REDUCE  FROM  RENVONT 
CA±L  REDUCE(PINELNSAVE4) 

PINELB  -  PINEL 

RETURN 

END 

SUBROUTINE  RIBAR(DASP) 

C  THIS  SUBROUTINE  MODIFIES  EXISTING  RIB  DIMENSIONS  TO  GIVE  A 
C  DESIRED  ASPECT  RATIO  WITHOUT  AFFECnNG  SHELL  INERTIA. 

C  H  «  RIB  HEIGHT,  W  »  RIB  WIDTH,  DASP  -  DESIRED  ASPECT  RATIO. 

COMMON  DUM(23).H,W 

COMMON  /COUNT/  ITERAT,INDIC,INDIC2,IND1C3, 

1  INDICS, NDESIGN, ORTHO, NOKG,NCURVE,OK,niERAM,OKM,nERAG,OKG,ISEC, 

2  NPLOT,NOKM,NRINGS,NEGLAM,INDIC7 
REAL  IE,IE1 

CALL  lECAUIE) 

C  PRINT  •,•  RIBARI  ■,DUM(8),DUM(9),DUM(1),DUM(.'’),H,W,DUM(5) 

DO  50  I  -  UOO 
ITERAT  -  ITERAT  +  1 
W  -  H/DASP 
CALL  lEC/UXIEI) 

Q  •  (IE-IE1)/AMAX1(IE,IEI) 

IF  (ABS(Q).LT..OOOI)  RETURN 
H  -  H/(L-Q/2.) 

50  CONTINUE 
CALL  IECAL(IE) 

RETURN 

END 

SUBROUTINE  RIBBIT(RINGS,SECDEP,DW,DASP) 

C  THIS  SUBROUTINE  WHIPS  UP  A  GOOD  INITIAL  DESIGN  OF  A  MONOCOOUE 
C  OR  RING  STIFFENED  SHELL  FROM  WHICH  THE  MANY  REFINEMENTS  ARE  MADE 
C  THROUGHOUT  THE  PROGRAM  IN  THE  DESIGN  OPTION. 

COMMON  AREA,IR,YR,LF,LS,ARM,LRO,T,GXY,EX,EY,EF,VXY,VYX,VPLAS,R, 

I  P,STR,SIGX,SIGY,DUMI(2),H,W,OD 
COMMON  /COUNT/  ITERAT,INDIC,INDIC2,INDIO, 
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1  INDIC5,NDESIGN,ORTHO.NOKG,NCURVE,OK,ITERAM,OICM,ITERAG,OKG,ISEC 

2  NPLOT,NOKM,NRJNGS,NEGLAM, INDIO 
INTEGER  ORTHO 

REAL  IR,K.L,LFa^ 

C  PRINT  •;  RIBBnr,RO.TAREA,YR,H.W.LS 

SIGMAY  -  SQRT(SIGX**2  +  SIOY'*2  -  SIGX'SIGY) 

ETX  -  EX 
ETY«  EY 
ETF  -  EF 
VLT-  VXY 
VTL-  VYX 

IF  (NCURVE.EQ.4)  GO  TO  10 

C  ASSUME80PERCENTYIELD  AND  ASTRESSRAT10  0F2. 

S  «  .8 
K  ■  J 

SF  -  S'SIGMAY 
SY-  SF 
SX  ■  K*SY 

IF  (ORTHO.EOjO)  THEN 
SX  -  SQRT(SX‘*2  +  SY**2  -  SX*SY) 

SY  -  SX 
END  IF 

C  PRINT  *,  •  CALL  MODULI  FROM  RIBBFT 

CALL  MODULI(SXSY.SF.ETX,ETY.ETF.VLT,VTL,NCURVE) 

10  CONTINUE 
RO  •  OD/2. 

NRINGS  -  RINGS 
LS  •  U(RINGS  +  1.) 

LF-  LS 

T  ■  P*RO/(2ySQRT(3.)*SIOMAY  +  P/1) 

DO  100  ITERAL  «  l.iOO 
C  PRINT  •,  •  RIBBITl  •.RO.R.TOF.VLT.VTL 

CALL  LOCAL<ETX,ETV’.VLT,VTL4»CRITUNCRrrL) 

C  PRINT  •,  •  RIBBIT2  ’.RO.R.T,  LFFCRITL,NCRITL,ETX,ETY 
IF(ITERAL.E0.1  AND.PCRrrLOE.P)  GO  TO  115 
F  -  (p.PCRrfL)/AMAXl(P,PCRITL) 

IF  (ABS(F).LT,.001)  GO  TO  113 
T  -  T/(l.  •  F/3.) 

100  CONTINUE 

113  IF  (NRINOS.OT.O)  GO  TO  116 
W-  0, 

H-  0. 
iAREA-  0. 

IR-  0. 

YR-  0. 

DW-  0. 

DASP-  0. 

SECDEP-0. 

RETURN 
116  CONTINUE 
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H  »  ABS(SECDEP  -  T) 

IF  (H.LT.T)  H  -  T 
W-  H/l 
NCRITL  •  -100 
YR  -  T+  H/2. 

ITERAG  -  0 

200  ITERAT  •  ITERAT  1 
AREA  -  WH 
IR  -  AREA*H**2/12. 

I  c  _  I  p  u/ 

ITERAG  •  ITERAO  +  1 
IF  (ITERAG  .EQ.  150)  GO  TO  207 
IF(W,LE.LS)  GO  TO  205 
LS  »  LF/2. 

CALL  RIBDW(LS) 

205  NCRITL  -  -100 

CALL  GENERL(ETX,ETV.ETF.VLT.VTL,PCRITG,NCRITL) 

G-  (P-PCRnG)/AMAXl(P4*CRrTG) 

IF  (ABS(G).LT..01)  GO  TO  207 
W-  W/(I.^) 

GO  TO  200 
207  CONTINUE 

IF  (DW.GT.0.  .AND.  W.OT.DW)  CALL  RIBDW(DW) 

IF  (ISECEQ.I)  CALL  RIBH(SECDEP) 

AR  -  H/W 

IF  (AR.LT.DASP)  CALL  RIBAR(DASP) 

AR  •  H/W 

IF  (AR.GT.4.AND.W.LT..125)  CALL  RIBDW(0.125) 

AREA  ■  W*H 
IR  -  AREA/I2.*H**2 
YR  -  T  +  H/2. 

LS  -  LF  -  W 

C  PRINT  •;  ribbita*,ro,tarea.yr,h,w,ls 
IF  (W.LE.LS)  RETURN 
LS  -  LF/2. 

CALL  RIBDW(LS) 

LS  -  LF  -  W 

C  PRINT  RIBBITB".RO,T,AREA,YR,H.W4-S 
RETURN 
END 

SUBROUTINE  RIBDW(DW) 

C  THIS  SUBROUTINE  MODIFIES  EXISTING  RIB  DIMENSIONS  TO  GIVE  A  DESIRED 
C  RIB  WIDTH  WITHOUT  /UTECTING  SHELL  INERTIA. 

C  H  -  RIB  HEIGHT,  W  -  RIB  WIDTH,  DW  -  DESIRED  WIDTH. 

COMMON  DUM(23),H,W 

COMMON  /COUNT/  ITERAT, INDIC,INDIC2,IND1C3, 

1  INDia,NDESIGN,ORTHO.NOKO.NCURVE,OK,ITERAM,01CM,ITERAG.OKG,ISEC, 

2  NPLOT,NOKM,NRINGS,NEGLAM,INDIC'7 
REAL  IE,IEI 

CALL  IECAL(IE) 
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C  PRINT  *,■  RIBDWI  •,DUM(8)4>UM(9)J)UM(I),DUM(3),H,W,DUM(5),DW 
W  -  DW 
DO  50  I  -  1300 
ITERAT  -  HERAT  +  1 
CALL  IECAL(IEI) 

Q  «  (1E1E1)/AMAXUIE.IE1) 

IF  (ABS(O).LT..0(X)l)  RETURN 
H  -  H/(1.-CW.) 

50  CONTINUE 
CALL  lECAMIE) 

RETURN 
END 

SUBROUTINE  RIBH(SECDEP) 

THIS  SUBROUTINE  MODIFIES  EXISTING  RIB  DIMENSIONS  TO  GIVE  A 
DESIRED  SHELL  SECHON  DEPTH  WITHOUT  AFFECTING  SHELL  INERTIA. 

T  -  SHELL  THICKNESS.  H  -  RIB  HEIGHT.  W  -  RIB  WIDTH. 

SECDEP  -  H  +  T  -  SECTION  DEPTH. 

COMMON  DUM  1(3),LF.LS,DUM2(3),TJ)UM3(I4).H.W 
COMMON  /COUNT/  ITERAT.INDIC,INDIC2.INDIC3, 

1  INDia.NDESION,ORTHO.NOKO.NCURVE.OK,ITERAM.OKM,ITERAG.OKG.ISEC. 

2  NPLOT.NOKM.NRINGS.NEGLAM.INDIC:7 
REAL  IE,IE1.LF.LS 
WW  «  W 
HH  •  H 

C  IF  (T.LEX).)  PRINT  *,  *  RIBH  •.H.W.U.LS.T.ITCRAT 
CALL  IECAL(IE) 

C  PRINT  RIBH  I  •.DUM2(3).T,DUM I(I).DUMI(3)JI.WXS 
H  -  ABS(SECDEP  -  T) 

DO  50  I  ■  1300 
ITERAT  -  ITERAT  I 
CALL  lECAUIEI) 

Q  -  (IE  IEI)/AMAX1(IE.IE1> 

C  IF(H.LT.O.)  PRINT  KlBHIE  •.H.W.LS.LF.IE.IEI.Q.I 
IF  (ABS(O).LT..0001)  RETURN 
W  -  W/(l..Q/2.) 

LS-  LF.  W 
IF  (LS.LT.W)  GO  TO  70 
50  (X)NnNUE 
CALL  IECAL(IE) 

RETURN 
70  CONTINUE 
H  -  HH 
m  WW 

CALL  IECAL(IE) 

LS-  LF/2. 

CALL  RIBDW(LS) 

LS  -  LF  .  W 

C  PRINT  •;  RIBH  B  •,DUM2(3),T,DUMI(I),DUM1(3),H.W.LS 
RETURN 
END 


119 


APPENDIX  B:  THESIS  CODE  INPUT  HLE 


The  program  input  file  allows  the  user  to  generally  describe  fixed  shell 
geometiy  variables  (such  as  L  and  OD),  properties  of  the  shell  and  ring  material(s), 
and  the  applied  load:  external  hydrostatic  pressure  (PSTR).  This  appendix  shows  the 
input  file  to  be  read  by  the  THESIS  code.  The  input  file  is  named  DAPS3  INP. 
This  is  the  same  input  file  read  by  the  DAPS3  code  when  used  as  a  "stand  alone" 
shell  design  program  (in  both  "design"  and  "analysis"  modes).  Adherence  to  the 
character  spacing  is  mandatory  to  insure  that  the  input  parameters  are  read  properly. 
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70  7$  60 


0  S 

1  1 

10  IS 

1  1 

20  25 

1  1 

80  38 

1  1 

(1) 

N1»  NE** 

N3"  N4' 

It 

(2) 

TITLE 

(«) 

HHH 

(4) 

All  A22 

011  022 

(5) 

EX 

VYX 

VPLA8 

SIOX 

(«) 

P8TR 

8I0Y 

810X1 

SIQYI 

(7) 

L 

00 

T 

RNO 

(•) 

AR 

IR 

YR 

LF 

(») 

RINGS 

8EC0EP 

ON 

0A8P 

(10) 

SIOXNON 

SlOYNON 

NX'*  NY'* 

NF'* 

(11)  6TM8X(1),  I  ■  I.NXPNT 

(12)  STRlNXd),  1  •  I.NXPNT 

(13)  6TR£8Y(1),  1  ••  1  ,NYFNT 

(14)  STRINYd).  1  ■  I.NYPMT 

(15)  8TnE8F(I),  1  •  I.NFPNr 

(16)  8T6INF(I),  I  •  I.NFFMT 


40  45  60  88  60  65 

I  I  I  I  I  I  I  I  I 


OXY 

EY 

EF 

PEL 

8TRX 

8TRY 

N8AVE 

RNOF 

U 

ARH 

RlNGSl 

RIN08F 

ORINQS 

I8EC 

INPUT  LtNC  FONMATS: 


<M 

|2 

i3> 

i4i 

iSi 

(6 

\l 

!  il 
01 


4tS 
7A10 
7A10 
0F10.S 
6F10.3 
7F10.8,  IS 
8F10.3 
OFIO.S 
7F10.3,  15 
:  2F10.3,  315 
THRU  (10):  10F8.0 


NOTE:  (1)  BLANK  ENTRIES  ARE  READ  A8  EERO 

(2  NOT  ALL  INPUT  LINES  REQUIRED,  OEPENOINO  ON  VALUES  OF  NPLOT,  NCURVE, 
N0E81QN,  AND  ORTHO 


•N1’  •  NPLOT 
*N2*  >  NCURVE 
*N3*  •  NDESION 
•N4*  •  ORTHO 
•NX’  •  NXPNT 
•NY*  •  NYPNT 
•NF*  •  NFPNT 
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APPENDIX  C:  THESIS  CODE  FLOW  CHART  LOGIC 


This  appendix  shows  the  flow  chart  logic  of  the  interaction  between  the 
THESIS  and  DAPS3  codes. 
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